As a user, you can use AstroFV to perform simulation of an astrophysical application. The following simulations are currently supported. Please click on the relevant one:
Introduction
Astrophysical scenarios involving relativistic flows occur in several phenomenon, most notably in the jets in extragalactic radio sources associated with active galactic nuclei. Similar to shock wave phenomenon in classical Newtonian fluid mechanics, strong shocks are a common feature in such astrophysical scenarios. Numerical simulations are therefore needed to study the formation, evolution and interaction of shock waves in relativistic fluids.
Relativistic shock tube problem (RSTP) involves the decay of an initially discontinuous two fluids into three elementary wave structures: a shock, a contact, and a rarefaction wave.
Governing Equations
This application considers the hydrodynamics model of a perfect fluid. The model is thus given by the hyperbolic system of conservation laws of the relativistic hydrodynamics.
Following Hujeirat and Font, the 1-dimensional hydrodynamical system in Minkowski space can be reformulated in terms of dynamical variables as:
where we define
= Relativistic Density
= Momentum
= Total energy (assuming conservation of mechanical energy)
= Transport Velocity
in Minkowski space
= four-velocity of the fluid
= rest mass density in locally intertial reference frame
= relativistic specific enthalpy
= specific internal energy
= pressure of ideal gas
The above equations are the desired equations which can be solved for the following cases:
-
Adiabatic conditions
- isothermal gas (specified with gamma = 0)
- non-isothermal gas with gamma = User defined (practically 4/3, 5/3, 2)
-
Initial conditions
= User input
= User input
= User input
= User input
= User input
-
Remarks
-
The hydrodynamical system represents a coupled set of non-linear advection equations
-
If fluid is considered isothermal, h~1 and the equation in E should be ignored (need not solved)
-
For non-isothermal case, Dh = D +
-
Dh ~ D for isothermal case
-
-
-
-
-
-
-
To get
we have used that
and
-
Numerical techniques
Using first order upwinding with constant grid size
Constant grid size is used
We use cell centered FV discretization for Density and Momentum conservation and vertex centered FV discretization for Total energy. This is shown below
For vertex centered approach, the transport velocity at FV cell surface is taken as the mean value across the boundary. i.e
The upwind flux across the FV cell surfaces is defined as:
The paramater controls the implicit scheme. Value 1.0 implies Implicit Euler. Other useful values like 0.5 which imply Crank-Nicholson are for future use.
The 3 set of equations can be written as:
which when re-arranged gives a tridiagonal system of equations as follows:
, then
and , then
- Remarks
- The dependence on values from (n+1)th iteration for Vx and P are handled by using an iterative procedure. This is explained in pseudo code below:
Expectation: We are in iteration n, given values of step n we want to calculate values for n+1 Input: Dn, Vn, Mn, Edn, Pn as the values from step n, iter_count as number of sub-iterations Procedure: Set: V* = Vn P* = Pn Loop: Execute the below sub-iteration for iter_count: #Solve the 3 equations D' = f1(Dn,Vn,V*) M' = f2(Mn,Vn,V*,P*) E' = f3(En,Vn,V*) #Update phase V* = update_Vx(D',M') P* = update_pressure(D',E') Repeat loop Finally: Vn+1 = V* Pn+1 = P* Dn+1 = D' Mn+1 = M' En+1 = E'
Key Data structures
- Parameters
class RSTPExplicitParams(Params):
def __init__(self,ncells,gamma=0,cfl=1.0):
Params.__init__(self,"explicit")
self.gamma = gamma
self.cfl = cfl
self.ncells = ncells
self.fv_boundary_strategy = None
class RSTPImplicitParams(Params):
def __init__(self,ncells,alpha,gamma=0,iter_count=1,cfl=0.50):
Params.__init__(self,"implicit")
self.gamma = gamma
self.ncells = ncells
self.alpha = alpha
self.fv_boundary_strategy = None
self.cfl = cfl #This is used to define time step size
self.iter_count = iter_countThese are used to define parameters for Explicit Euler time integration and Implicit time integration.
gamma = 0 defines isothermal gas
Time limits are default [0,1] and spatial limits are default [0,1].
Spatial grid size is automatically chosen by dividing spatial limit by ncells.
Temporal grid size is automatically derived from spatial grid size by using CFL parameter. delta_t = self.params.cfl*delta_x
Example:
eparams = RSTPExplicitParams(1000,0,0.25) #isothermal gas, CFL=0.25
iparams = RSTPImplicitParams(1000,1.0,0.0,2,0.7) #isothermal gas, Implicit Euler solver, iter_count=2, cfl=0.7 - Initial Values
class RSTPIV(InitialValues):
def __init__(self,Vx,Mx,D,Rho):
self.Vx = Vx
self.Mx = Mx
self.D = D
self.Rho = RhoUsed to provide initial values. The initial values are provided as a list for the two discontinuous fluids.
Example:
iv = RSTPIV(Vx=[0,0],Mx=[0,0],D=[1,10**-2],Rho=[1,10**-2])Usage examples
Number of FV cells = 1000, gamma = 4/3, CFL = 0.25
eparams = RSTPExplicitParams(1000,4/3,0.25)
eparams.set_fig_path('./figs/')
eparams.fv_boundary_strategy = FVTransverse #Default
eiv = RSTPIV(Vx=[0,0],Mx=[0,0],D=[1,10**-2],Rho=[1,10**-2])
ebv = RSTPBV()
test_explicit = RSTPTest(1,eparams,eiv,ebv,ode_strategy=ODEExplicit)
test_explicit.solve()Number of FV cells = 1000, gamma = 4/3, CFL = 0.7, no of subiterations = 2
iparams = RSTPImplicitParams(1000,1.0,4/3,2,0.7)
iparams.set_fig_path('./figs/')
iparams.fv_boundary_strategy = FVTransverse #Default
iiv = RSTPIV(Vx=[0,0],Mx=[0,0],D=[1,10**-2],Rho=[1,10**-2])
ibv = RSTPBV()
test_implicit = RSTPTest(2,iparams,iiv,ibv,ode_strategy=ODEImplicit)
test_implicit.solve()The output is generated as a set of 4 image files in the same directory as provided in set_fig_path() API.
The 4 image files correspond to plot of spatial variable x against Relativistic Density (D), Pressure(P), Lorentz factor () and Transport Velocity (
) respectively
An example is shown below:
