PLUTO
init.c File Reference

Collection of relativistic MHD shock-tube problems. More...

#include "pluto.h"
Include dependency graph for init.c:

Go to the source code of this file.

Functions

void Init (double *v, double x1, double x2, double x3)
 
void Analysis (const Data *d, Grid *grid)
 
void UserDefBoundary (const Data *d, RBox *box, int side, Grid *grid)
 

Detailed Description

Collection of relativistic MHD shock-tube problems.

This directory contains several 1D shock-tube configurations commonly used for benchmarking of numerical methods. There are 16 input parameters that can be used to specify entirely the left and right states with respect to the discontinuity initially placed at x=1/2. For the left state one has:

  • g_inputParam[RHO_LEFT]: density for the left state
  • g_inputParam[VX_LEFT]: normal (x) velocity for the left state
  • g_inputParam[VY_LEFT]: tangential (y) velocity for the left state
  • g_inputParam[VZ_LEFT]: tangential (z) velocity for the left state
  • g_inputParam[BY_LEFT]: tangential (y) magnetic field for the left state
  • g_inputParam[BZ_LEFT]: tangential (z) magnetic field for the left state
  • g_inputParam[PR_LEFT]: pressure in the left state

Similar parameters are used to specify the right states. In addition,

  • g_inputParam[BX_CONST]: normal component of magnetic field (this cannot have jump)
  • g_inputParam[GAMMA_EOS]: the specific heat ratio

The IDEAL equation of state is used. The available configurations are taken from

Conf Reference
01-021st problem in [MB06]
03-042nd problem in [MB06]
05-063rd problem in [MB06]
07-084th problem in [MB06]
09-102nd problem in [MUB09]
11-124th problem (generic Alfven test) in [MUB09]
13-14[MMB05]
rmhd_sod.11.jpg
Results for the generic Alfven test at t = 0.5 on 1600 grid zones for conf. #11.
Authors
A. Mignone (migno.nosp@m.ne@p.nosp@m.h.uni.nosp@m.to.i.nosp@m.t)
Date
Oct 6, 2014 16, 2014

References

  • [MB06]: "An HLLC Riemann solver for relativistic flows -- II. Magnetohydrodynamics", Mignone & Bodo MNRAS(2006) 368,1040
  • [MUB09]: "A five-wave Harten-Lax-van Leer Riemann solver for relativistic magnetohydrodynamics", Mignone et al., MNRAS(2009) 393,1141
  • [MMB05]: "Relativistic MHD simulations of jets with toroidal magnetic fields", Mignone, Massaglia & Bodo, SSRv (2005), 121,21
  • Balsara, ApJS (2001), 132, 83
  • Del Zanna, Bucciantini, Londrillo, A&A (2003)

Definition in file init.c.

Function Documentation

void Analysis ( const Data d,
Grid grid 
)

Perform runtime data analysis.

Parameters
[in]dthe PLUTO Data structure
[in]gridpointer to array of Grid structures

Compute volume-integrated magnetic pressure, Maxwell and Reynolds stresses. Save them to "averages.dat"

Definition at line 92 of file init.c.

97 {
98 
99 }
void Init ( double *  v,
double  x1,
double  x2,
double  x3 
)

The Init() function can be used to assign initial conditions as as a function of spatial position.

Parameters
[out]va pointer to a vector of primitive variables
[in]x1coordinate point in the 1st dimension
[in]x2coordinate point in the 2nd dimension
[in]x3coordinate point in the 3rdt dimension

The meaning of x1, x2 and x3 depends on the geometry:

\[ \begin{array}{cccl} x_1 & x_2 & x_3 & \mathrm{Geometry} \\ \noalign{\medskip} \hline x & y & z & \mathrm{Cartesian} \\ \noalign{\medskip} R & z & - & \mathrm{cylindrical} \\ \noalign{\medskip} R & \phi & z & \mathrm{polar} \\ \noalign{\medskip} r & \theta & \phi & \mathrm{spherical} \end{array} \]

Variable names are accessed by means of an index v[nv], where nv = RHO is density, nv = PRS is pressure, nv = (VX1, VX2, VX3) are the three components of velocity, and so forth.

Definition at line 62 of file init.c.

68 {
70 
71  if (x1 < 0.5){
73  v[VX1] = g_inputParam[VX_LEFT];
74  v[VX2] = g_inputParam[VY_LEFT];
75  v[VX3] = g_inputParam[VZ_LEFT];
77  v[BX2] = g_inputParam[BY_LEFT];
78  v[BX3] = g_inputParam[BZ_LEFT];
79  v[PRS] = g_inputParam[PR_LEFT];
80  }else{
84  v[VX3] = g_inputParam[VZ_RIGHT];
87  v[BX3] = g_inputParam[BZ_RIGHT];
88  v[PRS] = g_inputParam[PR_RIGHT];
89  }
90 }
#define BZ_RIGHT
double g_gamma
Definition: globals.h:112
#define PR_RIGHT
#define VX2
Definition: mod_defs.h:29
#define PR_LEFT
#define RHO
Definition: mod_defs.h:19
#define BY_LEFT
#define VY_RIGHT
#define VX1
Definition: mod_defs.h:28
#define RHO_LEFT
#define RHO_RIGHT
#define VZ_LEFT
#define VY_LEFT
#define BZ_LEFT
double g_inputParam[32]
Array containing the user-defined parameters.
Definition: globals.h:131
#define BX3
Definition: mod_defs.h:27
#define BX_CONST
#define VX_LEFT
#define VZ_RIGHT
#define GAMMA_EOS
#define BY_RIGHT
#define BX1
Definition: mod_defs.h:25
#define VX3
Definition: mod_defs.h:30
#define BX2
Definition: mod_defs.h:26
#define VX_RIGHT
void UserDefBoundary ( const Data d,
RBox box,
int  side,
Grid grid 
)

Assign user-defined boundary conditions.

Parameters
[in,out]dpointer to the PLUTO data structure containing cell-centered primitive quantities (d->Vc) and staggered magnetic fields (d->Vs, when used) to be filled.
[in]boxpointer to a RBox structure containing the lower and upper indices of the ghost zone-centers/nodes or edges at which data values should be assigned.
[in]sidespecifies the boundary side where ghost zones need to be filled. It can assume the following pre-definite values: X1_BEG, X1_END, X2_BEG, X2_END, X3_BEG, X3_END. The special value side == 0 is used to control a region inside the computational domain.
[in]gridpointer to an array of Grid structures.

Assign user-defined boundary conditions in the lower boundary ghost zones. The profile is top-hat:

\[ V_{ij} = \left\{\begin{array}{ll} V_{\rm jet} & \quad\mathrm{for}\quad r_i < 1 \\ \noalign{\medskip} \mathrm{Reflect}(V) & \quad\mathrm{otherwise} \end{array}\right. \]

where $ V_{\rm jet} = (\rho,v,p)_{\rm jet} = (1,M,1/\Gamma)$ and M is the flow Mach number (the unit velocity is the jet sound speed, so $ v = M$).

Assign user-defined boundary conditions:

  • left side (x beg): constant shocked values
  • bottom side (y beg): constant shocked values for x < 1/6 and reflective boundary otherwise.
  • top side (y end): time-dependent boundary: for $ x < x_s(t) = 10 t/\sin\alpha + 1/6 + 1.0/\tan\alpha $ we use fixed (post-shock) values. Unperturbed values otherwise.

Assign user-defined boundary conditions at inner and outer radial boundaries. Reflective conditions are applied except for the azimuthal velocity which is fixed.

Assign user-defined boundary conditions.

Parameters
[in/out]d pointer to the PLUTO data structure containing cell-centered primitive quantities (d->Vc) and staggered magnetic fields (d->Vs, when used) to be filled.
[in]boxpointer to a RBox structure containing the lower and upper indices of the ghost zone-centers/nodes or edges at which data values should be assigned.
[in]sidespecifies on which side boundary conditions need to be assigned. side can assume the following pre-definite values: X1_BEG, X1_END, X2_BEG, X2_END, X3_BEG, X3_END. The special value side == 0 is used to control a region inside the computational domain.
[in]gridpointer to an array of Grid structures.

Set the injection boundary condition at the lower z-boundary (X2-beg must be set to userdef in pluto.ini). For $ R \le 1 $ we set constant input values (given by the GetJetValues() function while for $ R > 1 $ the solution has equatorial symmetry with respect to the z=0 plane. To avoid numerical problems with a "top-hat" discontinuous jump, we smoothly merge the inlet and reflected value using a profile function Profile().

Provide inner radial boundary condition in polar geometry. Zero gradient is prescribed on density, pressure and magnetic field. For the velocity, zero gradient is imposed on v/r (v = vr, vphi).

The user-defined boundary is used to impose stress-free boundary and purely vertical magnetic field at the top and bottom boundaries, as done in Bodo et al. (2012). In addition, constant temperature and hydrostatic balance are imposed. For instance, at the bottom boundary, one has:

\[ \left\{\begin{array}{lcl} \DS \frac{dp}{dz} &=& \rho g_z \\ \noalign{\medskip} p &=& \rho c_s^2 \end{array}\right. \qquad\Longrightarrow\qquad \frac{p_{k+1}-p_k}{\Delta z} = \frac{p_{k} + p_{k+1}}{2c_s^2}g_z \]

where $g_z$ is the value of gravity at the lower boundary. Solving for $p_k$ at the bottom boundary where $k=k_b-1$ gives:

\[ \left\{\begin{array}{lcl} p_k &=& \DS p_{k+1} \frac{1-a}{1+a} \\ \noalign{\medskip} \rho_k &=& \DS \frac{p_k}{c_s^2} \end{array}\right. \qquad\mathrm{where}\qquad a = \frac{\Delta z g_z}{2c_s^2} > 0 \]

where, for simplicity, we keep constant temperature in the ghost zones rather than at the boundary interface (this seems to give a more stable behavior and avoids negative densities). A similar treatment holds at the top boundary.

Assign user-defined boundary conditions. At the inner boundary we use outflow conditions, except for velocity which we reset to zero when there's an inflow

Definition at line 101 of file init.c.

106 {
107 }