File paralleltransform.hxx#
-
class ParallelTransform#
- #include <paralleltransform.hxx>
Calculates the values of a field along the magnetic field (y direction)
This is used inside Mesh to represent a variety of possible numerical schemes, including Shifted Metric and FCI
Subclassed by FCITransform, ParallelTransformIdentity, ShiftedMetric, ShiftedMetricInterp
Public Functions
-
virtual ~ParallelTransform() = default#
-
virtual void calcParallelSlices(Field3D &f) = 0#
Given a 3D field, calculate and set the Y up down fields.
-
inline virtual void integrateParallelSlices(Field3D &f)#
Calculate Yup and Ydown fields by integrating over mapped points This should be used for parallel divergence operators
-
virtual Field3D toFieldAligned(const Field3D &f, const std::string ®ion = "RGN_ALL") = 0#
Convert a field into field-aligned coordinates so that the y index is along the magnetic field
-
virtual Field3D fromFieldAligned(const Field3D &f, const std::string ®ion = "RGN_ALL") = 0#
Convert back from field-aligned coordinates into standard form
-
inline virtual Field2D toFieldAligned(const Field2D &f, const std::string ®ion = "RGN_ALL")#
Field2D are axisymmetric, so transformation to or from field-aligned coordinates is a null operation.
-
virtual bool canToFromFieldAligned() const = 0#
-
inline virtual std::vector<PositionsAndWeights> getWeightsForYUpApproximation(int i, int j, int k)#
-
inline virtual std::vector<PositionsAndWeights> getWeightsForYDownApproximation(int i, int j, int k)#
-
inline virtual std::vector<PositionsAndWeights> getWeightsForYApproximation(int i, int j, int k, int yoffset)#
-
inline virtual void outputVars([[maybe_unused]] Options &output_options)#
Output variables used by a ParallelTransform instance to
output_options
.
-
virtual bool requiresTwistShift(bool twist_shift_enabled, YDirectionType ytype) = 0#
If
twist_shift_enabled
is true, does aField3D
with Y directionytype
require a twist-shift at branch cuts on closed field lines?
Protected Functions
-
virtual void checkInputGrid() = 0#
This method should be called in the constructor to check that if the grid has a ‘parallel_transform’ variable, it has the correct value
Protected Attributes
-
Options &options#
Options for this ParallelTransform.
-
struct PositionsAndWeights#
-
virtual ~ParallelTransform() = default#
-
class ParallelTransformIdentity : public ParallelTransform#
- #include <paralleltransform.hxx>
This class implements the simplest form of ParallelTransform where the domain is a logically rectangular domain, and yup() and ydown() refer to the same field.
Public Functions
-
virtual void calcParallelSlices(Field3D &f) override#
Merges the yup and ydown() fields of f, so that f.yup() = f.ydown() = f
-
inline virtual Field3D toFieldAligned(const Field3D &f, const std::string ®ion = "RGN_ALL") override#
The field is already aligned in Y, so this does nothing
-
inline virtual FieldPerp toFieldAligned(const FieldPerp &f, const std::string ®ion = "RGN_ALL") override#
-
inline virtual Field3D fromFieldAligned(const Field3D &f, const std::string ®ion = "RGN_ALL") override#
The field is already aligned in Y, so this does nothing
-
inline virtual FieldPerp fromFieldAligned(const FieldPerp &f, const std::string ®ion = "RGN_ALL") override#
-
inline virtual std::vector<PositionsAndWeights> getWeightsForYApproximation(int i, int j, int k, int yoffset) override#
-
inline virtual bool canToFromFieldAligned() const override#
-
inline virtual bool requiresTwistShift(bool twist_shift_enabled, YDirectionType ytype) override#
If
twist_shift_enabled
is true, does aField3D
with Y directionytype
require a twist-shift at branch cuts on closed field lines?
Protected Functions
-
virtual void checkInputGrid() override#
This method should be called in the constructor to check that if the grid has a ‘parallel_transform’ variable, it has the correct value
-
virtual void calcParallelSlices(Field3D &f) override#
-
class ShiftedMetric : public ParallelTransform#
- #include <paralleltransform.hxx>
Shifted metric method Each Y location is shifted in Z with respect to its neighbours so that the grid is orthogonal in X-Z, but requires interpolation to calculate the values of points along field-lines.
In this implementation the interpolation is done using FFTs in Z
Public Functions
-
ShiftedMetric() = delete#
-
ShiftedMetric(Mesh &mesh, CELL_LOC location, Field2D zShift, BoutReal zlength_in, Options *opt = nullptr)#
-
virtual void calcParallelSlices(Field3D &f) override#
Calculates the yup() and ydown() fields of f by taking FFTs in Z and applying a phase shift.
-
virtual Field3D toFieldAligned(const Field3D &f, const std::string ®ion = "RGN_ALL") override#
Uses FFTs and a phase shift to align the grid points with the y coordinate (along magnetic field usually).
Note that the returned field will no longer be orthogonal in X-Z, and the metric tensor will need to be changed if X derivatives are used.
Shift the field so that X-Z is not orthogonal, and Y is then field aligned.
-
virtual FieldPerp toFieldAligned(const FieldPerp &f, const std::string ®ion = "RGN_ALL") override#
-
virtual Field3D fromFieldAligned(const Field3D &f, const std::string ®ion = "RGN_ALL") override#
Converts a field back to X-Z orthogonal coordinates from field aligned coordinates.
Shift back, so that X-Z is orthogonal, but Y is not field aligned.
-
virtual FieldPerp fromFieldAligned(const FieldPerp &f, const std::string ®ion = "RGN_ALL") override#
-
inline virtual std::vector<PositionsAndWeights> getWeightsForYApproximation(int i, int j, int k, int yoffset) override#
-
inline virtual bool canToFromFieldAligned() const override#
-
inline virtual bool requiresTwistShift(bool twist_shift_enabled, YDirectionType ytype) override#
If
twist_shift_enabled
is true, does aField3D
with Y directionytype
require a twist-shift at branch cuts on closed field lines?
Protected Functions
-
virtual void checkInputGrid() override#
This method should be called in the constructor to check that if the grid has a ‘parallel_transform’ variable, it has the correct value
Private Functions
-
inline Field2D shiftZ(const Field2D &f, const Field2D &zangle, const std::string region = "RGN_NOX") const#
Shift a 2D field in Z. Since 2D fields are constant in Z, this has no effect
-
Field3D shiftZ(const Field3D &f, const Tensor<dcomplex> &phs, const YDirectionType y_direction_out, const std::string ®ion = "RGN_NOX") const#
Shift a 3D field or FieldPerp
f
by the given phasephs
in ZCalculates FFT in Z, multiplies by the complex phase and inverse FFTS.
- Parameters:
f – [in] The field to shift
phs – [in] The phase to shift by
y_direction_out – [in] The value to set yDirectionType of the result to
-
FieldPerp shiftZ(const FieldPerp &f, const Tensor<dcomplex> &phs, const YDirectionType y_direction_out, const std::string ®ion = "RGN_NOX") const#
-
void shiftZ(const BoutReal *in, const dcomplex *phs, BoutReal *out) const#
Shift a given 1D array, assumed to be in Z, by the given
zangle
- Parameters:
in – [in] A 1D array of length mesh.LocalNz
phs – [in] Phase shift, assumed to have length (mesh.LocalNz/2 + 1) i.e. the number of modes
out – [out] A 1D array of length mesh.LocalNz, already allocated
-
void cachePhases()#
Calculate and store the phases for to/from field aligned and for the parallel slices using zShift
-
std::vector<Field3D> shiftZ(const Field3D &f, const std::vector<ParallelSlicePhase> &phases) const#
Shift a 3D field
f
in Z to all the parallel slices inphases
- Parameters:
f – [in] The field to shift
phases – [in] The phase and offset information for each parallel slice
- Returns:
The shifted parallel slices
Private Members
-
CELL_LOC location = {CELL_CENTRE}#
-
Field2D zShift#
This is the shift in toroidal angle (z) which takes a point from X-Z orthogonal to field-aligned along Y.
-
int nmodes#
-
Tensor<dcomplex> toAlignedPhs#
Cache of phase shifts for transforming from X-Z orthogonal coordinates to field-aligned coordinates
-
Tensor<dcomplex> fromAlignedPhs#
Cache of phase shifts for transforming from field-aligned coordinates to X-Z orthogonal
-
std::vector<ParallelSlicePhase> parallel_slice_phases#
Cache of phase shifts for the parallel slices. Slices are stored in the following order: {+1, …, +n, -1, …, -n} slice[i] stores offset i+1 slice[n + i] stores offset -(i+1) where i goes from 0 to (n-1), with n the number of y guard cells
-
struct ParallelSlicePhase#
coordinates
Helper POD for parallel slice phase shifts
-
ShiftedMetric() = delete#