PLUTO
init.c File Reference

Sedov-Taylor blast wave. More...

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

Go to the source code of this file.

Functions

void Init (double *us, 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

Sedov-Taylor blast wave.

Set the initial condition for a Sedov-Taylor blast wave problem. Ambient density and pressure are constant everywhere and equal to rho0 and p0. While rho0 is an input parameter, p0=1.e-5. A large amount of energy is deposited in a small volume (a line, circle or sphere depending on geometry and dimensions) consisting of 2 or 3 computational zones.

The input parameters read from pluto.ini are labeled as:

  • g_inputParam[ENRG0]: initial energy of the blast wave (not energy density);
  • g_inputParam[DNST0]: initial constant density;
  • g_inputParam[GAMMA]: specific heat ratio.

The available configurations refer to:

  1. Cartesian (1D)
  2. Cylindrical (1D)
  3. Spherical (1D)
  4. Cartesian (2D, equivalent to cylindrical)
  5. Cartesian (3D, equivalent to spherical)
hd_sedov.02.jpg
Density plot at t=0.5 for configuration #02 (cylindrical blast)
Author
A. Mignone (migno.nosp@m.ne@p.nosp@m.h.uni.nosp@m.to.i.nosp@m.t)
Date
July 08, 2014

References:

  • L.I. Sedov, "Similarity and Dimensional Methods in Mechanics", Academic Press, New York, 1959.

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

Definition at line 94 of file init.c.

99 {
100 
101 }
void Init ( double *  us,
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 40 of file init.c.

46 {
47  double dr, vol, r;
48 
50 
51 /* --------------------------------------------------
52  dr is the size of the initial energy deposition
53  region: 2 ghost zones.
54  -------------------------------------------------- */
55 
56  #if DIMENSIONS == 1
57  dr = 2.0*(g_domEnd[IDIR]-g_domBeg[IDIR])/(double)NX1;
58  #else
59  dr = 3.5*(g_domEnd[IDIR]-g_domBeg[IDIR])/(double)NX1;
60  #endif
61 
62 /* ---------------------------------------
63  compute region volume
64  --------------------------------------- */
65 
66  #if (GEOMETRY == CARTESIAN) && (DIMENSIONS == 1)
67  vol = 2.0*dr;
68  #elif (GEOMETRY == CYLINDRICAL && DIMENSIONS == 1)|| \
69  (GEOMETRY == CARTESIAN && DIMENSIONS == 2)
70  vol = CONST_PI*dr*dr;
71  #elif (GEOMETRY == SPHERICAL && DIMENSIONS == 1)|| \
72  (GEOMETRY == CYLINDRICAL && DIMENSIONS == 2)|| \
73  (GEOMETRY == CARTESIAN && DIMENSIONS == 3)
74  vol = 4.0/3.0*CONST_PI*dr*dr*dr;
75  #else
76  print1 ("! Init: geometrical configuration not allowed\n");
77  QUIT_PLUTO(1);
78  #endif
79 
80  r = EXPAND(x1*x1, + x2*x2, +x3*x3);
81  r = sqrt(r);
82 
83  us[RHO] = g_inputParam[DNST0];
84  us[VX1] = 0.0;
85  us[VX2] = 0.0;
86  us[VX3] = 0.0;
87 
88  if (r <= dr) us[PRS] = (g_gamma - 1.0)*g_inputParam[ENRG0]/vol;
89  else us[PRS] = 1.e-5;
90 
91  us[TRC] = 0.0;
92 }
double g_gamma
Definition: globals.h:112
#define VX2
Definition: mod_defs.h:29
long int NX1
Number of interior zones in the X1 directions (boundaries excluded) for the local processor...
Definition: globals.h:48
void print1(const char *fmt,...)
Definition: amrPluto.cpp:511
#define ENRG0
#define RHO
Definition: mod_defs.h:19
#define DNST0
#define TRC
Definition: pluto.h:581
#define VX1
Definition: mod_defs.h:28
#define GAMMA
#define IDIR
Definition: pluto.h:193
double g_inputParam[32]
Array containing the user-defined parameters.
Definition: globals.h:131
double g_domBeg[3]
Lower limits of the computational domain.
Definition: globals.h:125
#define VX3
Definition: mod_defs.h:30
#define CONST_PI
.
Definition: pluto.h:269
double g_domEnd[3]
Upper limits of the computational domain.
Definition: globals.h:126
#define QUIT_PLUTO(e_code)
Definition: macros.h:125

Here is the call graph for this function:

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.

Definition at line 104 of file init.c.

109 { }