PLUTO  4.4-patch2
Functions
eigenv.c File Reference

Wave-speeds and characteristic decomposition for the MHD equations. More...

#include "pluto.h"

Functions

void MaxSignalSpeed (const State *state, double *cmin, double *cmax, int beg, int end)
 
void Eigenvalues (double **v, double *csound2, double **lambda, int beg, int end)
 
void PrimEigenvectors (const State *state, int beg, int end)
 
void ConsEigenvectors (double *u, double *v, double a2, double **Lc, double **Rc, double *lambda)
 
void PrimToChar (double **Lp, double *v, double *w)
 

Detailed Description

This file contains various functions containing Jacobian-related information such as characteristic signal speeds, eigenvalues and eigenvector decomposition for the MHD module.

The function MaxSignalSpeed() computes the maximum and minimum characteristic signal velocity for the MHD equations.

The function Eigenvalues() computes the 7 characteristic waves.

The function PrimEigenvectors() calculates left and right eigenvectors and the corresponding eigenvalues for the primitive form the the MHD equations with an adiabatic or isothermal EoS.

The function ConsEigenvectors() provides the characteristic decomposition of the convervative MHD equations.

The function PrimToChar() compute the matrix-vector multiplcation between the L matrix (containing primitive left eigenvectors) and the vector v. The result containing the characteristic variables is stored in the vector w = L.v

Authors
A. Mignone (migno.nosp@m.ne@t.nosp@m.o.inf.nosp@m.n.it)
P. Tzeferacos
Date
Oct 24, 2019

Function Documentation

◆ ConsEigenvectors()

void ConsEigenvectors ( double *  u,
double *  v,
double  a2,
double **  Lc,
double **  Rc,
double *  lambda 
)

Provide conservative eigenvectors for MHD equations.

Parameters
[in]uarray of conservative variables
[in]varray of primitive variables
[in]a2square of sound speed
[out]Lcleft conservative eigenvectors
[out]Rcright conservative eigenvectors
[out]lambdaeigenvalues

Reference

  • "High-order conservative finite difference GLM-MHD schemes for cell-centered MHD"
    Mignone, Tzeferacos & Bodo, JCP (2010) 229, 5896
  • "A High-order WENO Finite Difference Scheme for the Equations of Ideal MHD"
    Jiang,Wu, JCP 150,561 (1999)

With the following corrections:

The components (By, Bz) of L_{1,7} page 571 should be +{rho}a instead of -. Also, since Bx is not considered as a parameter, one must also add a component in the fast and slow waves. This can be seen by forming the conservative eigenvectors from the primitive ones, see the paper from Powell, Roe Linde.

The Lc_{1,3,5,7}[Bx] components, according to Serna 09 are -(1-g_gamma)*a_{f,s}*Bx and not (1-g_gamma)*a_{f,s}*Bx. Both are orthonormal though. She is WRONG! -Petros-

◆ Eigenvalues()

void Eigenvalues ( double **  v,
double *  csound2,
double **  lambda,
int  beg,
int  end 
)

Compute eigenvalues for the MHD equations

Parameters
[in]v1-D array of primitive variables
[out]csound21-D array containing the square of sound speed
[out]lambda1-D array [i][nv] containing the eigenvalues
[in]begstarting index of computation
[in]endfinal index of computation

◆ MaxSignalSpeed()

void MaxSignalSpeed ( const State state,
double *  cmin,
double *  cmax,
int  beg,
int  end 
)

Compute the maximum and minimum characteristic velocities for the MHD equations, cmin= v - cf, cmax = v + cf

Parameters
[in]statepointer to a State structure
[out]cmin1-D array containing the leftmost characteristic
[out]cmin1-D array containing the rightmost characteristic
[in]begstarting index of computation
[in]endfinal index of computation

◆ PrimEigenvectors()

void PrimEigenvectors ( const State state,
int  beg,
int  end 
)

Provide left and right eigenvectors and corresponding eigenvalues for the PRIMITIVE form of the MHD equations (adiabatic & isothermal cases).

Parameters
[in,out]statePointer to State structure
[in]begstarting grid index
[in]endfinal grid index
Note
Eigenvectors state->Lp and state->Rp must be initialized to zero before since only non-zero entries are treated here.

Wave names and their order are defined as enumeration constants in mod_defs.h. Notice that, the characteristic decomposition may differ depending on the way div.B is treated.

Advection modes associated with passive scalars are simple cases for which lambda = u (entropy mode) and l = r = (0, ... , 1, 0, ...). For this reason they are NOT defined here and must be treated seperately.

References:

  • "Notes on the Eigensystem of Magnetohydrodynamics",
    P.L. Roe, D.S. Balsara, SIAM Journal on Applied Mathematics, 56, 57 (1996)
  • "A solution adaptive upwind scheme for ideal MHD",
    K. G. Powell, P. L. Roe, and T. J. Linde, Journal of Computational Physics, 154, 284-309, (1999).
  • "A second-order unsplit Godunov scheme for cell-centered MHD: the CTU-GLM scheme"
    Mignone & Tzeferacos, JCP (2010) 229, 2117
  • "ATHENA: A new code for astrophysical MHD", J. Stone, T. Gardiner, ApJS, 178, 137 (2008)

The derivation of the isothermal eigenvectors follows the consideration given in roe.c

◆ PrimToChar()

void PrimToChar ( double **  Lp,
double *  v,
double *  w 
)

Compute the matrix-vector multiplcation between the the L matrix (containing primitive left eigenvectors) and the vector v. The result containing the characteristic variables is stored in the vector w = L.v

For efficiency purpose, multiplication is done explicitly, so that only nonzero entries of the left primitive eigenvectors are considered.

Parameters
[in]LpLeft eigenvectors
[in]v(difference of) primitive variables
[out]w(difference of) characteristic variables