File fv_ops.hxx#

namespace FV

Functions

inline Field3D Div_a_Laplace_perp(const Field3D &a, const Field3D &f)#
template<typename CellEdges = MC>
Field3D Div_par(const Field3D &f_in, const Field3D &v_in, const Field3D &wave_speed_in, bool fixflux = true)#

Finite volume parallel divergence

Preserves the sum of f*J*dx*dy*dz over the domain

NB: Uses to/from FieldAligned coordinates

Parameters:
  • f_in[in] The field being advected. This will be reconstructed at cell faces using the given CellEdges method

  • v_in[in] The advection velocity. This will be interpolated to cell boundaries using linear interpolation

  • wave_speed_in[in] Local maximum speed of all waves in the system at each

  • fixflux[in] Fix the flux at the boundary to be the value at the midpoint (for boundary conditions)

template<typename CellEdges = MC>
Field3D Div_f_v(const Field3D &n_in, const Vector3D &v, bool bndry_flux)#

Div ( n * v ) &#8212; Magnetic drifts

This uses the expression

Div( A ) = 1/J * d/di ( J * A^i )

Hence the input vector should be contravariant

Note: Uses to/from FieldAligned

template<typename CellEdges = MC>
Field3D Div_par_mod(const Field3D &f_in, const Field3D &v_in, const Field3D &wave_speed_in, Field3D &flow_ylow, bool fixflux = true)#

Finite volume parallel divergence

NOTE: Modified version, applies limiter to velocity and field Performs better (smaller overshoots) than Div_par

Preserves the sum of f*J*dx*dy*dz over the domain

Parameters:
  • f_in[in] The field being advected. This will be reconstructed at cell faces using the given CellEdges method

  • v_in[in] The advection velocity. This will be interpolated to cell boundaries using linear interpolation

  • wave_speed_in[in] Local maximum speed of all waves in the system at each

  • fixflux[in] Fix the flux at the boundary to be the value at the midpoint (for boundary conditions)

  • flow_ylow[out] Flow at the lower Y cell boundary Already includes area factor * flux

template<typename CellEdges = MC>
Field3D Div_par_fvv(const Field3D &f_in, const Field3D &v_in, const Field3D &wave_speed_in, bool fixflux = true)#

This operator calculates Div_par(f v v) It is used primarily (only?) in the parallel momentum equation.

This operator is used rather than Div(f fv) so that the values of f and v are consistent with other advection equations: The product fv is not interpolated to cell boundaries.

struct Stencil1D#
#include <fv_ops.hxx>

Stencil used for Finite Volume calculations which includes cell face values L and R

Public Members

BoutReal c#

Cell centre values.

BoutReal m#
BoutReal p#
BoutReal mm = BoutNaN#
BoutReal pp = BoutNaN#
BoutReal L = BoutNaN#

Left cell face value.

BoutReal R = BoutNaN#

Right cell face value.

struct Upwind#
#include <fv_ops.hxx>

First order upwind for testing

Public Functions

inline void operator()(Stencil1D &n)#
struct Fromm#
#include <fv_ops.hxx>

Fromm method

Public Functions

inline void operator()(Stencil1D &n)#
struct MinMod#
#include <fv_ops.hxx>

Second order slope limiter method

Limits slope to minimum absolute value of left and right gradients. If at a maximum or minimum slope set to zero, i.e. reverts to first order upwinding

Public Functions

inline void operator()(Stencil1D &n)#

Private Static Functions

static inline BoutReal _minmod(BoutReal a, BoutReal b)#

Internal helper function for minmod slope limiter

If the inputs have different signs then returns zero, otherwise chooses the value with the minimum magnitude.

struct MC#
#include <fv_ops.hxx>

Monotonised Central (MC) second order slope limiter (Van Leer)

Limits the slope based on taking the slope with the minimum absolute value from central, 2*left and 2*right. If any of these slopes have different signs then the slope reverts to zero (i.e. 1st-order upwinding).

Public Functions

inline void operator()(Stencil1D &n)#

Private Static Functions

static inline BoutReal minmod(BoutReal a, BoutReal b, BoutReal c)#
struct Superbee#
#include <fv_ops.hxx>

Superbee limiter

This corresponds to the limiter function φ(r) = max(0, min(2r, 1), min(r,2)

The value at cell right (i.e. i + 1/2) is:

n.R = n.c - φ(r) (n.c - (n.p + n.c)/2) = n.c + φ(r) (n.p - n.c)/2

Four regimes: a) r < 1/2 -> φ(r) = 2r n.R = n.c + gL b) 1/2 < r < 1 -> φ(r) = 1 n.R = n.c + gR/2 c) 1 < r < 2 -> φ(r) = r n.R = n.c + gL/2 d) 2 < r -> φ(r) = 2 n.R = n.c + gR

where the left and right gradients are: gL = n.c - n.m gR = n.p - n.c

Public Functions

inline void operator()(Stencil1D &n)#