PLUTO
Functions
tvdlf.c File Reference

Lax-Friedrechs (Rusanov) Riemann solver for HD. More...

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

Go to the source code of this file.

Functions

void LF_Solver (const State_1D *state, int beg, int end, double *cmax, Grid *grid)
 

Detailed Description

Lax-Friedrechs (Rusanov) Riemann solver for HD.

Solve the Riemann problem for the adiabatic and isothermal HD equations using the Lax-Friedrichs Rusanov Riemann solver with local maximum characteristic speed:

\[ \hat{F}_{i+\HALF} = \frac{F^L_{i+\HALF} + F^R_{i+\HALF}}{2} - c^{\rm max}_{i+\HALF}\frac{U^R_{i+\HALF} - U^L_{i-\HALF}}{2} \qquad{\rm where}\quad c^{\rm max} = \lambda_{\rm max}\left( \frac{V^R_{i+\HALF} + V^L_{i+\HALF}}{2}\right) \]

where $\lambda_{\rm \max}$ is a function of the arithmetic average between the left and the right states.

On input, this function takes left and right primitive state vectors state->vL and state->vR at zone edge i+1/2; On output, return flux and pressure vectors at the same interface i+1/2 (note that the i refers to i+1/2).

Also during this step, compute maximum wave propagation speed (cmax) for explicit time step computation.

Reference:

  • "Comparison of Some Flux Correced Transport and Total Variation Diminishing Numerical Schemes for Hydrodynamics and Magnetohydrodynamics Problems", Toth and Odstrcil, JCP (1996), 128,82
Authors
A. Mignone (migno.nosp@m.ne@p.nosp@m.h.uni.nosp@m.to.i.nosp@m.t)
Date
April 7, 2014

Definition in file tvdlf.c.

Function Documentation

void LF_Solver ( const State_1D *  state,
int  beg,
int  end,
double *  cmax,
Grid *  grid 
)
Parameters
[in,out]statepointer to State_1D structure
[in]beginitial grid index
[out]endfinal grid index
[out]cmax1D array of maximum characteristic speeds
[in]gridpointer to array of Grid structures.

Definition at line 38 of file tvdlf.c.

49 {
50  int nv, i;
51  static double **fL, **fR, **vRL;
52  static double *cRL_min, *cRL_max, *pL, *pR, *a2L, *a2R;
53  double *vR, *vL, *uR, *uL, *flux;
54 
55  if (fR == NULL){
56  fR = ARRAY_2D(NMAX_POINT, NFLX, double);
57  fL = ARRAY_2D(NMAX_POINT, NFLX, double);
58  vRL = ARRAY_2D(NMAX_POINT, NFLX, double);
59  cRL_min = ARRAY_1D(NMAX_POINT, double);
60  cRL_max = ARRAY_1D(NMAX_POINT, double);
61  pR = ARRAY_1D(NMAX_POINT, double);
62  pL = ARRAY_1D(NMAX_POINT, double);
63  a2R = ARRAY_1D(NMAX_POINT, double);
64  a2L = ARRAY_1D(NMAX_POINT, double);
65  }
66 
67 /* ----------------------------------------------------
68  compute sound speed & fluxes at zone interfaces
69  ---------------------------------------------------- */
70 
71  SoundSpeed2 (state->vL, a2L, NULL, beg, end, FACE_CENTER, grid);
72  SoundSpeed2 (state->vR, a2R, NULL, beg, end, FACE_CENTER, grid);
73 
74  Flux (state->uL, state->vL, a2L, fL, pL, beg, end);
75  Flux (state->uR, state->vR, a2R, fR, pR, beg, end);
76 
77 /* ---------------------------------------------------------------------
78  Compute average state in order to get the local max characteristic
79  velocity.
80  In order to avoid underestimating this speed at strong reflecting
81  boundaries (or when vxL = -vxR in general) , we average the
82  absolute value of the normal velocity.
83  ------------------------------------------------------------------- */
84 
85  for (i = beg; i <= end; i++){
86  VAR_LOOP(nv) vRL[i][nv] = 0.5*(state->vL[i][nv] + state->vR[i][nv]);
87  vRL[i][VXn] = 0.5*(fabs(state->vL[i][VXn]) + fabs(state->vR[i][VXn]));
88  }
89 
90 /* ---------------------------------------------------------------------
91  Compute sound speed, max and min signal speeds
92  --------------------------------------------------------------------- */
93 
94  SoundSpeed2 (vRL, a2R, NULL, beg, end, FACE_CENTER, grid);
95  MaxSignalSpeed (vRL, a2R, cRL_min, cRL_max, beg, end);
96 
97  for (i = beg; i <= end; i++) {
98  uR = state->uR[i];
99  uL = state->uL[i];
100  flux = state->flux[i];
101 
102  /* -- compute max eigenvalue -- */
103 
104  cmax[i] = MAX(fabs(cRL_max[i]), fabs(cRL_min[i]));
105  g_maxMach = MAX(g_maxMach, fabs(vRL[i][VXn])/sqrt(a2R[i]));
106 
107  /* -- compute fluxes -- */
108 
109  for (nv = NFLX; nv--; ) {
110  flux[nv] = 0.5*(fL[i][nv] + fR[i][nv] - cmax[i]*(uR[nv] - uL[nv]));
111  }
112  state->press[i] = 0.5*(pL[i] + pR[i]);
113  }
114 }
MAX
#define MAX(a, b)
Definition: macros.h:101
STATE_1D::flux
double ** flux
upwind flux computed with the Riemann solver
Definition: structs.h:149
Flux
void Flux(double **u, double **w, double *a2, double **fx, double *p, int beg, int end)
Definition: fluxes.c:23
STATE_1D::vR
double ** vR
Primitive variables to the right of the interface, .
Definition: structs.h:139
SoundSpeed2
void SoundSpeed2(double **v, double *cs2, double *h, int beg, int end, int pos, Grid *grid)
Definition: eos.c:16
g_maxMach
double g_maxMach
The maximum Mach number computed during integration.
Definition: globals.h:119
VAR_LOOP
#define VAR_LOOP(n)
Definition: macros.h:226
NFLX
#define NFLX
Definition: mod_defs.h:32
FACE_CENTER
#define FACE_CENTER
Definition: pluto.h:206
VXn
int VXn
Definition: globals.h:73
STATE_1D::uR
double ** uR
same as vR, in conservative vars
Definition: structs.h:145
ARRAY_1D
#define ARRAY_1D(nx, type)
Definition: prototypes.h:170
NMAX_POINT
long int NMAX_POINT
Maximum number of points among the three directions, boundaries excluded.
Definition: globals.h:62
i
int i
Definition: analysis.c:2
STATE_1D::vL
double ** vL
Primitive variables to the left of the interface, .
Definition: structs.h:136
STATE_1D::press
double * press
Upwind pressure term computed with the Riemann solver.
Definition: structs.h:164
MaxSignalSpeed
void MaxSignalSpeed(double **v, double *cs2, double *cmin, double *cmax, int beg, int end)
Definition: eigenv.c:34
ARRAY_2D
#define ARRAY_2D(nx, ny, type)
Definition: prototypes.h:171
STATE_1D::uL
double ** uL
same as vL, in conservative vars
Definition: structs.h:144

Here is the call graph for this function: