PLUTO  4.0
 All Data Structures Files Functions Variables Enumerations Macros Pages
Enumerations | Functions
mod_defs.h File Reference

Set labels, indexes and prototypes for the MHD module. More...

Go to the source code of this file.

Macros

Vector Potential Labels

These may only be used in the STARTUP / INIT functions. They're convenient in obtaining a discretization that preserve the divergence-free condition (for staggered field) or if you simply wish to initialize the magnetic field from the vector potential.

#define AX1   (NVAR + 1)
 
#define AX2   (NVAR + 2)
 
#define AX3   (NVAR + 3)
 

Enumerations

enum  KWAVES
 

Functions

void BackgroundField (double x1, double x2, double x3, double *B0)
 
int ConsToPrim (double **, real **, int, int, unsigned char *)
 
void Eigenvalues (double **, double *, double **, int, int)
 
void PrimEigenvectors (double *, double, double, double *, double **, double **)
 
void ConsEigenvectors (double *, double *, double, double **, double **, double *)
 
void Enthalpy (double **, double *, int, int)
 
void Entropy (double **, double *, int, int)
 
void Flux (double **, double **, double *, double **, double **, double *, int, int)
 
void MaxSignalSpeed (double **, double *, double *, double *, double **, int, int)
 
void PrimToCons (double **, double **, int, int)
 
void PrimRHS (double *, double *, double, double, double *)
 
void PrimSource (const State_1D *, int, int, double *, double *, double **, Grid *)
 

Detailed Description

Contains basic macro definitions, structure definitions and global variable declarations used by the MHD module.

Author
A. Mignone (migno.nosp@m.ne@p.nosp@m.h.uni.nosp@m.to.i.nosp@m.t)
Date
Sep, 10, 2012

Enumeration Type Documentation

enum KWAVES

Label the different waves in increasing order following the number of vector components.

IMPORTANT: the KPSI_GLMM & KPSI_GLMP modes are present only in the MHD-GLM formulation. We keep them at the END of the enumeration so we can skip them in unnecessary loops. Please do NOT change them !

Function Documentation

void BackgroundField ( double  x1,
double  x2,
double  x3,
double *  B0 
)

Define the component of a static, curl-free background magnetic field.

Parameters
[in]x1position in the 1st coordinate direction $x_1$
[in]x2position in the 2nd coordinate direction $x_2$
[in]x3position in the 3rd coordinate direction $x_3$
[out]B0array containing the vector componens of the background magnetic field
void ConsEigenvectors ( double *  u,
double *  v,
double  a2,
double **  Lc,
double **  Rc,
double *  lambda 
)

Provide conservative eigenvectors for HD equations.

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

References:

  • "Riemann Solvers and Numerical Methods for Fluid Dynamics",
    Toro, 1997 Springer-Verlag (Page 107)

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-

int ConsToPrim ( double **  ucons,
double **  uprim,
int  ibeg,
int  iend,
unsigned char *  flag 
)

Convert from conservative to primitive variables.

Parameters
[in]uconsarray of conservative variables
[out]uprimarray of primitive variables
[in]begstarting index of computation
[in]endfinal index of computation
[out]flagarray of flags tagging zones where conversion went wrong.
Returns
Return (0) if conversion was succesful in every zone [ibeg,iend]. Otherwise, return a non-zero integer number giving the bit flag(s) turned on during the conversion process. In this case, flag contains the failure codes of those zones where where conversion did not go through.

Convert from conservative to primitive variables.

Parameters
[in]uconsarray of conservative variables
[out]uprimarray of primitive variables
[in]begstarting index of computation
[in]endfinal index of computation
[out]flagarray of flags tagging zones where conversion went wrong.
Returns
This function returns 0 on success. Otherwise, a non-zero integer number giving the bit flag turned on during the conversion.

Convert from conservative to primitive variables.

Parameters
[in]uconsarray of conservative variables
[out]uprimarray of primitive variables
[in]begstarting index of computation
[in]endfinal index of computation
[out]flagarray of flags tagging zones where conversion went wrong.
Returns
Return (0) if conversion was succesful in every zone [ibeg,iend]. Otherwise, return a non-zero integer number giving the bit flag(s) turned on during the conversion process. In this case, flag contains the failure codes of those zones where where conversion did not go through.

Convert from conservative to primitive variables.

Parameters
[in]uconsarray of conservative variables
[out]uprimarray of primitive variables
[in]begstarting index of computation
[in]endfinal index of computation
[out]flagarray of flags tagging zones where conversion went wrong.
Returns
This function returns 0 on success. Otherwise, a non-zero integer number giving the bit flag turned on during the conversion.
void Eigenvalues ( double **  v,
double *  csound2,
double **  lambda,
int  beg,
int  end 
)

Compute eigenvalues for the HD 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

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
void Enthalpy ( double **  v,
double *  h,
int  beg,
int  end 
)

Compute the enthalpy.

Parameters
[in]v1D array of primitive quantities
[in]h1D array of enthalpy values
[in]beginitial index of computation
[in]endfinal index of computation
Returns
This function has no return value.
void Entropy ( double **  v,
double *  s,
int  beg,
int  end 
)

Compute the entropy.

Parameters
[in]v1D array of primitive quantities
[in]s1D array of entropy values
[in]isinitial index of computation
[in]iefinal index of computation
Returns
This function has no return value.

Compute the entropy.

Parameters
[in]v1D array of primitive quantities
[in]s1D array of entropy values
[in]beginitial index of computation
[in]endfinal index of computation
Returns
This function has no return value.
void Flux ( double **  ucons,
double **  wprim,
double *  a2,
double **  bck,
double **  fx,
double *  p,
int  beg,
int  end 
)
Parameters
[in]ucons1D array of conserved quantities
[in]wprim1D array of primitive quantities
[in]a21D array of sound speeds
[in]bck1D array of background field values
[out]fx1D array of fluxes (total pressure excluded)
[out]p1D array of pressure values
[in]beginitial index of computation
[in]endfinal index of computation
Returns
This function has no return value.
void MaxSignalSpeed ( double **  v,
double *  cs2,
double *  cmin,
double *  cmax,
double **  bgf,
int  beg,
int  end 
)

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

Parameters
[in]v1-D array of primitive variables
[out]cs21-D array containing the square of the sound speed
[out]cmin1-D array containing the leftmost characteristic
[out]cmin1-D array containing the rightmost characteristic
[in]bgf1-D array containing the background magnetic field
[in]begstarting index of computation
[in]endfinal index of computation
void PrimEigenvectors ( double *  q,
double  a2,
double  h,
double *  lambda,
double **  LL,
double **  RR 
)

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

Parameters
[in]qvector of primitive variables
[in]a2sound speed squared
[in]henthalpy
[out]lambdaeigenvalues
[out]LLleft primitive eigenvectors
[out]RRright primitive eigenvectors
Note
It is highly recommended that LL and RR 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.

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:

  • "Riemann Solvers and Numerical Methods for Fluid Dynamics",
    Toro, 1997 Springer-Verlag

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

Parameters
[in]qvector of primitive variables
[in]a2sound speed squared
[in]henthalpy
[out]lambdaeigenvalues
[out]LLleft primitive eigenvectors
[out]RRright primitive eigenvectors
Note
It is highly recommended that LL and RR 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

void PrimRHS ( real *  w,
real *  dw,
real  cs2,
real  h,
real *  Adw 
)

Compute the matrix-vector multiplication $ A(\mathbf{v})\cdot d\mathbf{v} $ where A is the matrix of the quasi-linear form of the HD equations.

Parameters
[in]wvector of primitive variables;
[in]dwlimited (linear) slopes;
[in]cs2local sound speed;
[in]hlocal enthalpy;
[out]Adwmatrix-vector product.
Returns
This function has no return value.

Compute the matrix-vector multiplication $ A(\mathbf{v})\cdot d\mathbf{v} $ where A is the matrix of the quasi-linear form of the MHD equations.

References

  • "A solution adaptive upwind scheme for ideal MHD", Powell et al., JCP (1999) 154, 284
  • "An unsplit Godunov method for ideal MHD via constrained transport" Gardiner & Stone, JCP (2005) 205, 509
Parameters
[in]vvector of primitive variables
[in]dvlimited (linear) slopes
[in]cs2local sound speed
[in]hlocal enthalpy
[out]AdVmatrix-vector product
Returns
This function has no return value.

Compute the matrix-vector multiplication $ A(\mathbf{v})\cdot d\mathbf{v} $ where A is the matrix of the quasi-linear form of the HD equations.

Parameters
[in]wvector of primitive variables;
[in]dwlimited (linear) slopes;
[in]cs2local sound speed;
[in]hlocal enthalpy;
[out]Adwmatrix-vector product.
Returns
This function has no return value.

Compute the matrix-vector multiplication $ A(\mathbf{v})\cdot d\mathbf{v} $ where A is the matrix of the quasi-linear form of the MHD equations.

References

  • "A solution adaptive upwind scheme for ideal MHD", Powell et al., JCP (1999) 154, 284
  • "An unsplit Godunov method for ideal MHD via constrained transport" Gardiner & Stone, JCP (2005) 205, 509
Parameters
[in]vvector of primitive variables
[in]dvlimited (linear) slopes
[in]cs2local sound speed
[in]hlocal enthalpy
[out]AdVmatrix-vector product
Returns
This function has no return value.
Note
In the 7-wave and 8-wave formulations we use the the same matrix being decomposed into right and left eigenvectors during the Characteristic Tracing step. Note, however, that it DOES NOT include two additional terms (-vy*dV[BX] for By, -vz*dv[BX] for Bz) that are needed in the

7-wave form and are added using source terms.

Note
In the 7-wave and 8-wave formulations we use the the same matrix being decomposed into right and left eigenvectors during the Characteristic Tracing step. Note, however, that it DOES NOT include two additional terms (-vy*dV[BX] for By, -vz*dv[BX] for Bz) that are needed in the

7-wave form and are added using source terms.

Note
In the 7-wave and 8-wave formulations we use the the same matrix being decomposed into right and left eigenvectors during the Characteristic Tracing step. Note, however, that it DOES NOT include two additional terms (-vy*dV[BX] for By, -vz*dv[BX] for Bz) that are needed in the

7-wave form and are added using source terms.

void PrimSource ( const State_1D state,
int  beg,
int  end,
double *  a2,
double *  h,
double **  src,
Grid grid 
)

Compute source terms of the HD equations in primitive variables.

  • Geometrical sources;
  • Gravity;
  • Fargo source terms.

The rationale for choosing during which sweep a particular source term has to be incorporated should match the same criterion used during the conservative update. For instance, in polar or cylindrical coordinates, curvilinear source terms are included during the radial sweep only.

Parameters
[in]statepointer to a State_1D structure;
[in]beginitial index of computation;
[in]endfinal index of computation;
[in]a2array of sound speed;
[in]harray of enthalpies (not needed in MHD);
[out]srcarray of source terms;
[in]gridpointer to a Grid structure.
Returns
This function has no return value.
Note
This function does not work in spherical coordinates yet.

Compute source terms of the MHD equations in primitive variables. These include:

  • Geometrical sources;
  • Gravity;
  • terms related to divergence of B control (Powell eight wave and GLM);
  • FARGO source terms.

The rationale for choosing during which sweep a particular source term has to be incorporated should match the same criterion used during the conservative update. For instance, in polar or cylindrical coordinates, curvilinear source terms are included during the radial sweep only.

Parameters
[in]statepointer to a State_1D structure
[in]beginitial index of computation
[in]endfinal index of computation
[in]a2array of sound speed
[in]harray of enthalpies (not needed in MHD)
[out]srcarray of source terms
[in]gridpointer to a Grid structure
Note
This function does not work in spherical coordinates yet. For future implementations we annotate hereafter the induction equation in spherical coordinates:

\[ \partial_tB_r + \frac{1}{r}\partial_\theta E_\phi - \frac{1}{r\sin\theta}\partial_\phi E_\theta = -E_\phi\cot\theta/r \]

\[ \partial_t B_\theta + \frac{1}{r\sin\theta}\partial_\phi E_r - \partial_rE_\phi = E_\phi/r \]

\[ \partial_t B_\phi + \partial_r E_\theta - \frac{1}{r}\partial_\theta E_r = - E_\theta/r\]

where

\[ E_\phi = -(v \times B)_\phi = - (v_r B_\theta - v_\theta B_r) \,,\qquad E_\theta = -(v \times B)_\theta = - (v_\phi B_r - v_r B_\phi) \]

Compute source terms of the HD equations in primitive variables.

  • Geometrical sources;
  • Gravity;
  • Fargo source terms.

The rationale for choosing during which sweep a particular source term has to be incorporated should match the same criterion used during the conservative update. For instance, in polar or cylindrical coordinates, curvilinear source terms are included during the radial sweep only.

Parameters
[in]statepointer to a State_1D structure;
[in]beginitial index of computation;
[in]endfinal index of computation;
[in]a2array of sound speed;
[in]harray of enthalpies (not needed in MHD);
[out]srcarray of source terms;
[in]gridpointer to a Grid structure.
Returns
This function has no return value.
Note
This function does not work in spherical coordinates yet.

Compute source terms of the MHD equations in primitive variables. These include:

  • Geometrical sources;
  • Gravity;
  • terms related to divergence of B control (Powell eight wave and GLM);
  • FARGO source terms.

The rationale for choosing during which sweep a particular source term has to be incorporated should match the same criterion used during the conservative update. For instance, in polar or cylindrical coordinates, curvilinear source terms are included during the radial sweep only.

Parameters
[in]statepointer to a State_1D structure
[in]beginitial index of computation
[in]endfinal index of computation
[in]a2array of sound speed
[in]harray of enthalpies (not needed in MHD)
[out]srcarray of source terms
[in]gridpointer to a Grid structure
Note
This function does not work in spherical coordinates yet. For future implementations we annotate hereafter the induction equation in spherical coordinates:

\[ \partial_tB_r + \frac{1}{r}\partial_\theta E_\phi - \frac{1}{r\sin\theta}\partial_\phi E_\theta = -E_\phi\cot\theta/r \]

\[ \partial_t B_\theta + \frac{1}{r\sin\theta}\partial_\phi E_r - \partial_rE_\phi = E_\phi/r \]

\[ \partial_t B_\phi + \partial_r E_\theta - \frac{1}{r}\partial_\theta E_r = - E_\theta/r\]

where

\[ E_\phi = -(v \times B)_\phi = - (v_r B_\theta - v_\theta B_r) \,,\qquad E_\theta = -(v \times B)_\theta = - (v_\phi B_r - v_r B_\phi) \]

void PrimToCons ( double **  uprim,
double **  ucons,
int  ibeg,
int  iend 
)

Convert primitive variables to conservative variables.

Parameters
[in]uprimarray of primitive variables
[out]uconsarray of conservative variables
[in]begstarting index of computation
[in]endfinal index of computation

Convert primitive variables in conservative variables.

Parameters
[in]uprimarray of primitive variables
[out]uconsarray of conservative variables
[in]begstarting index of computation
[in]endfinal index of computation

Convert primitive variables to conservative variables.

Parameters
[in]uprimarray of primitive variables
[out]uconsarray of conservative variables
[in]begstarting index of computation
[in]endfinal index of computation

Convert primitive variables in conservative variables.

Parameters
[in]uprimarray of primitive variables
[out]uconsarray of conservative variables
[in]begstarting index of computation
[in]endfinal index of computation