Richtmyer Meshkov Instability - mlwong/HAMeRS GitHub Wiki

Tutorial 05: Richtmyer-Meshkov Instability

Tutorial Goal

In this tutorial, the setup of inviscid and viscous simulations of Richtmyer-Meshkov 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
  • Multiple levels of adaptive mesh refinement using more than one type of sensors (value and gradient sensors) to identify regions for refinement
  • Sixth-order conservative diffusive/viscous flux method
  • Different models for mass diffusivity and viscosity

Problem Description

Richtmyer-Meshkov instability (RMI) is a hydrodynamic instability that appears when a shock wave passes through a material interface between fluids of different densities. The instability is initiated from the misalignment in the pressure and density gradients across the shock wave and material interface. This tutorial is a single-mode Richtmyer-Meshkov instability problem with domain size [0, 16λ] x [0, λ). The initial conditions of the problem is given by the following table:

Pre-shock air Post-shock air Pre-shock SF6
ρ (kg/m3) 1.145598 1.616874 5.972856
p (Pa) 101325 164859 101325
T (K) 298 343.53 298
u (m/s) 436.201332 309.060123 436.201332
v (m/s) 0 0 0

The initial perturbed interface is located at

The initial material interface has a smooth profile and is given by:

where v are any primitive variables near the initial interface and erf() is the error function. Subscripts L and R denote the left and right interface conditions. ΔD is the distance from the initial perturbed material interface. This problem is periodic in the vertical direction. Both left and right boundary conditions are extrapolated from interior solutions. We have chosen λ=1mm and εi=0.006/128mm.

The figure below shows the initial configuration of the problem:

Initial Conditions

There are two separate initial condition files for the inviscid and viscous versions of this problem. The file containing the initial conditions of the inviscid version of this problem (RichtmyerMeshkovInstability2D.cpp) can be found in problems/Euler/initial_conditions and the initial conditions should be set up by ln -sf <absolute path to RichtmyerMeshkovInstability2D.cpp> src/apps/Euler/EulerInitialConditions.cpp. The file containing the initial conditions of the viscous version of this problem (RichtmyerMeshkovInstability2D.cpp) can be found in problems/Navier-Stokes/initial_conditions.

Input File Configurations

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.

Application = "Euler"

or

Application = "Navier-Stokes"

Next, either the block Euler or Navier-Stokes is set depending on whether the problem is inviscid or viscous. Here we only consider the viscous case:

NavierStokes
{
    // Name of project
    project_name = "2D Richtmyer-Meshkov 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.09312, 1.39909
            species_R = 56.927, 296.803
        }

        // Equation of mass diffusivity to use
        equation_of_mass_diffusivity = "REID"

        Equation_of_mass_diffusivity_mixing_rules
        {
            species_epsilon_by_k = 222.1, 78.6
            species_sigma        = 5.128, 3.711
            species_M            = 146.055, 28.0135
        }

        // Equation of shear viscosity to use
        equation_of_shear_viscosity = "CHAPMAN_ENSKOG"

        Equation_of_shear_viscosity_mixing_rules
        {
            species_epsilon_by_k = 222.1, 78.6
            species_sigma        = 5.128, 3.711
            species_M            = 146.055, 28.0135
        }

        // Equation of bulk viscosity to use
        equation_of_bulk_viscosity = "CRAMER"

        Equation_of_bulk_viscosity_mixing_rules
        {
            species_gamma = 1.09312, 1.39909

            species_A_r = 0.0, -3.15e-5
            species_B_r = 0.0, 1.58e-7

            species_c_v_v = 7.73937, 0.0
            species_A_v   = 0.2064e-5, 0.0
            species_B_v   = 121.0, 0.0
            species_C_v   = -339.0, 0.0

            species_M = 146.055, 28.0135
        }

        // Equation of thermal conductivity to use
        equation_of_thermal_conductivity = "PRANDTL"

        Equation_of_thermal_conductivity_mixing_rules
        {
            species_c_p          = 668.286, 1040.50
            species_Pr           = 0.79, 0.71
            species_epsilon_by_k = 222.1, 78.6
            species_sigma        = 5.128, 3.711
            species_M            = 146.055, 28.0135

            equation_of_shear_viscosity = "CHAPMAN_ENSKOG"
        }
    }

    // Convective flux reconstructor to use
    convective_flux_reconstructor = "WCNS6_LD_HLLC_HLL"

    Convective_flux_reconstructor{}

    use_conservative_form_diffusive_flux = TRUE

    // Diffusive flux reconstructor to use
    diffusive_flux_reconstructor = "SIXTH_ORDER"

    Diffusive_flux_reconstructor{}

    Boundary_data
    {
        boundary_edge_xlo
        {
            boundary_condition = "FLOW"
        }
        boundary_edge_xhi
        {
            boundary_condition = "FLOW"
        }

        boundary_node_xlo_ylo
        {
            boundary_condition = "XFLOW"
        }
        boundary_node_xhi_ylo
        {
            boundary_condition = "XFLOW"
        }
        boundary_node_xlo_yhi
        {
            boundary_condition = "XFLOW"
        }
        boundary_node_xhi_yhi
        {
            boundary_condition = "XFLOW"
        }
    }

    Value_tagger
    {
        variables = "MASS_FRACTION"

        uses_global_tol_up = FALSE
        uses_global_tol_lo = FALSE

        uses_local_tol_up = TRUE
        uses_local_tol_lo = TRUE

        local_tol_up = 0.999
        local_tol_lo = 0.001
    }

    Gradient_tagger
    {
        gradient_sensors = "JAMESON_GRADIENT"

        JAMESON_GRADIENT
        {
            Jameson_gradient_variables = "PRESSURE"

            Jameson_gradient_tol = 0.002
        }
    }
}

Since the conservative four-equation model for multi-species flow is used, flow_model is set to "FOUR_EQN_CONSERVATIVE". One can also try to use constant mass diffusivities and viscosities by setting the mixing rules in Flow_model as:

        ...

        // Equation of mass diffusivity to use
        equation_of_mass_diffusivity = "CONSTANT"

        Equation_of_mass_diffusivity_mixing_rules
        {
            species_D = 9.76e-6, 9.76e-6
        }

        // Equation of shear viscosity to use
        equation_of_shear_viscosity = "CONSTANT"

        Equation_of_shear_viscosity_mixing_rules
        {
            species_mu = 1.53e-5, 1.81e-5
            species_M  = 146.055, 28.0135
        }

        // Equation of bulk viscosity to use
        equation_of_bulk_viscosity = "CONSTANT"

        Equation_of_bulk_viscosity_mixing_rules
        {
            species_mu_v = 0.0, 0.0
            species_M    = 146.055, 28.0135
        }

        // Equation of thermal conductivity to use
        equation_of_thermal_conductivity = "PRANDTL"

        Equation_of_thermal_conductivity_mixing_rules
        {
            species_c_p = 668.286, 1040.50
            species_Pr  = 0.79, 0.71
            species_mu  = 1.53-e5, 1.81e-5
            species_M   = 146.055, 28.0135

            equation_of_shear_viscosity = "CONSTANT"
        }

        ...

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), (511, 31)]
    x_lo         = 0.0, 0.0     // Lower end of computational domain
    x_up         = 0.016, 0.001 // Upper end of computational domain

    // Periodic dimension. A non-zero value indicates that the direction is periodic
    periodic_dimension = 0, 1
}

The following settings are used for ExtendedTagAndInitialize to use both value and gradient sensor to tag cells for refinement:

ExtendedTagAndInitialize
{
    at_0
    {
        cycle = 0
        // First tagging method for refinement to use
        tag_0
        {
            tagging_method = "VALUE_DETECTOR"
        }
        // Second tagging method for refinement to use
        tag_1
        {
            tagging_method = "GRADIENT_DETECTOR"
        }
    }
}

The parameters of adaptive mesh refinement are set in PatcHierarchy. In this tutorial, the maximum number of mesh levels is set to 3 and all levels have the same refinement ratio to coarser levels.

PatchHierarchy
{
    // Maximum number of levels in hierarchy
    max_levels = 3

    ratio_to_coarser
    {
        // Vector ratio to next coarser level
        level_1 = 2, 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             = 4.0e-5 // Final simulation time
    grow_dt              = 1.0e0  // Growth factor for timesteps
    max_integrator_steps = 10000  // Max number of simulation timesteps
    tag_buffer           = 2, 2   // array of integer values (one for each level that
                                  // may be refined representing the number of cells
                                  // by which tagged cells are buffered)
}

Running HAMeRS

To run the simulation with one core, first put the input file inside a directory named tutorial_05 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>

Visualization Using VisIt

The data can be visualized by opening dumps.visit inside the viz folders with VisIt.

The density field of inviscid simulation at t=4e-5s:

The density field of viscous simulation with physical diffusivity/viscosity models at t=4e-5s:

The density field of viscous simulation with constant diffusivity/viscosity at t=4e-5s:

⚠️ **GitHub.com Fallback** ⚠️