PLUTO  4.0
 All Data Structures Files Functions Variables Enumerations Macros Pages
Functions
prim_eqn.c File Reference

Compute the right hand side of the primitive MHD equations. More...

#include "pluto.h"

Functions

void PrimRHS (double *v, double *dv, double cs2, double h, double *Adv)
 
void PrimSource (const State_1D *state, int beg, int end, double *a2, double *h, double **src, Grid *grid)
 

Detailed Description

Implements the right hand side of the quasi-linear form of the MHD equations. In 1D this may be written as

\[ \partial_t{\mathbf{V}} = - A\cdot\partial_x\mathbf{V} + \mathbf{S} \]

where $ A $ is the matrix of the primitive form of the equations, $ S $ is the source term.

Reference

The function PrimRHS() implements the first term while PrimSource() implements the source term part.

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

Function Documentation

void PrimRHS ( double *  v,
double *  dv,
double  cs2,
double  h,
double *  Adv 
)

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.

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

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) \]