File field3d.hxx#

Defines

BOUT_FIELD3D_H#

Functions

FieldPerp operator+(const Field3D &lhs, const FieldPerp &rhs)#
FieldPerp operator-(const Field3D &lhs, const FieldPerp &rhs)#
FieldPerp operator*(const Field3D &lhs, const FieldPerp &rhs)#
FieldPerp operator/(const Field3D &lhs, const FieldPerp &rhs)#
Field3D operator+(const Field3D &lhs, const Field3D &rhs)#
Field3D operator-(const Field3D &lhs, const Field3D &rhs)#
Field3D operator*(const Field3D &lhs, const Field3D &rhs)#
Field3D operator/(const Field3D &lhs, const Field3D &rhs)#
Field3D operator+(const Field3D &lhs, const Field2D &rhs)#
Field3D operator-(const Field3D &lhs, const Field2D &rhs)#
Field3D operator*(const Field3D &lhs, const Field2D &rhs)#
Field3D operator/(const Field3D &lhs, const Field2D &rhs)#
Field3D operator+(const Field3D &lhs, BoutReal rhs)#
Field3D operator-(const Field3D &lhs, BoutReal rhs)#
Field3D operator*(const Field3D &lhs, BoutReal rhs)#
Field3D operator/(const Field3D &lhs, BoutReal rhs)#
Field3D operator+(BoutReal lhs, const Field3D &rhs)#
Field3D operator-(BoutReal lhs, const Field3D &rhs)#
Field3D operator*(BoutReal lhs, const Field3D &rhs)#
Field3D operator/(BoutReal lhs, const Field3D &rhs)#
Field3D operator-(const Field3D &f)

Unary minus. Returns the negative of given field, iterates over whole domain including guard/boundary cells.

Field3D pow(const Field3D &lhs, const Field2D &rhs, const std::string &rgn = "RGN_ALL")

Exponent: pow(lhs, lhs) is lhs raised to the power of rhs

Extra overloads not provided by the templates in field.hxx

This loops over the entire domain, including guard/boundary cells by default (can be changed using the rgn argument). If CHECK >= 3 then the result will be checked for non-finite numbers

FieldPerp pow(const Field3D &lhs, const FieldPerp &rhs, const std::string &rgn = "RGN_ALL")
void checkData(const Field3D &f, const std::string &region = "RGN_NOBNDRY")

Throw an exception if f is not allocated or if any elements are non-finite (for CHECK > 2). Loops over all points including the boundaries by default (can be changed using the rgn argument

Field3D filter(const Field3D &var, int N0, const std::string &rgn = "RGN_ALL")

Fourier filtering, removes all except one mode

Parameters:
  • var[in] Variable to apply filter to

  • N0[in] The component to keep

  • rgn[in] The region to calculate the result over

Field3D lowPass(const Field3D &var, int zmax, bool keep_zonal, const std::string &rgn = "RGN_ALL")

Fourier low pass filtering. Removes modes higher than zmax and optionally the zonal component

Parameters:
  • var[in] Variable to apply filter to

  • zmax[in] Maximum mode in Z

  • keep_zonal[in] Keep the zonal component if true

  • rgn[in] The region to calculate the result over

inline Field3D lowPass(const Field3D &var, int zmax, int keep_zonal, REGION rgn = RGN_ALL)#

The argument keep_zonal used to be integer “zmin” — this was a misnomer. Please use the version above which uses a bool instead

inline Field3D lowPass(const Field3D &var, int zmax, const std::string rgn = "RGN_ALL")#

Fourier low pass filtering. Removes modes higher than zmax

Parameters:
  • var[in] Variable to apply filter to

  • zmax[in] Maximum mode in Z

  • rgn[in] The region to calculate the result over

void shiftZ(Field3D &var, int jx, int jy, double zangle)

Perform a shift by a given angle in Z

Parameters:
  • var[inout] The variable to be modified in-place

  • jx[in] X index

  • jy[in] Y index

  • zangle[in] The Z angle to apply

void shiftZ(Field3D &var, BoutReal zangle, const std::string &rgn = "RGN_ALL")#

Apply a phase shift by a given angle zangle in Z to all points

Parameters:
  • var[inout] The variable to modify in-place

  • zangle[in] The angle to shift by in Z

  • rgn[in] The region to calculate the result over

Field2D DC(const Field3D &f, const std::string &rgn = "RGN_ALL")

Average in the Z direction

Parameters:
  • f[in] Variable to average

  • rgn[in] The region to calculate the result over

void invalidateGuards(Field3D &var)

Force guard cells of passed field var to NaN.

inline BOUT_HOST_DEVICE Field3D & ddt (Field3D &f)

Returns a reference to the time-derivative of a field f

Wrapper around member function f.timeDeriv()

template<>
inline std::string toString(const Field3D &val)#

toString template specialisation Defined in utils.hxx

bool operator==(const Field3D &a, const Field3D &b)

Test if two fields are the same, by calculating the minimum absolute difference between them

std::ostream &operator<<(std::ostream &out, const Field3D &value)

Output a string describing a Field3D to a stream.

class Field3D : public Field#
#include <field3d.hxx>

Class for 3D X-Y-Z scalar fields.

This class represents a scalar field defined over the mesh. It handles memory management, and provides overloaded operators for operations on the data, iterators and access methods.

Initialisation

Fields can be declared in any scope (even global), but cannot be accessed by index or used until the data is allocated.

Field3D f;   // Declare variable, no data allocated
f(0,0,0) = 1.0; // Error !

f = 0.0;  // Allocates memory, fills with value (0.0)

Field3D g(1.0); // Declares, allocates memory, fills with value (1.0)

Field3D h;   // not allocated
h.allocate();  // Data array allocated, values undefined
f(0,0,0) = 1.0; // ok
Copy-on-Write

A field is a reference to the underlying data array, so setting one field equal to another has the effect of making both fields share the same underlying data

Field3D f(0.0);
Field3D g = f; // f and g now share data
f(0,0,0) = 1.0; // g is also modified
Setting the entire field equal to a new value changes the reference:
Field3D f(0.0);
Field3D g = f; // f and g now share data
g = 1.0;   // g and f are now separate
To ensure that a field is unique, call allocate() which will make a copy of the underlying data if it is shared.
Field3D f(0.0);
Field3D g = f; // f and g now share data
g.allocate();  // Data copied so g and f don't share data
f(0,0,0) = 1.0; // ok
Data access

Individual data indices can be accessed by index using round brackets:

Field3D f;
f(0,1,2) = 1.0;  // Set value
BoutReal val = f(2,1,3);  // Get value
If CHECK is greater than 2, this function will perform bounds checking. This will significantly slow calculations.

Some methods, such as FFT routines, need access to a pointer to memory. For the Z dimension this can be done by passing only the X and Y indices

BoutReal *data = f(0,1);
data now points to f(0,1,0) and can be incremented to move in Z.

Iteration

To loop over all points in a field, a for loop can be used to get the indices:

Field3D f(0.0); // Allocate, set to zero

for( const auto &i : f ) {  // Loop over all points, with index i
  f[i] = 1.0;
}
There is also more explicit looping over regions:
for( const auto &i : f.region(RGN_ALL) ) {  // Loop over all points, with index i
  f[i] = 1.0;
}
Parallel (y) derivatives

In several numerical schemes the mapping along magnetic fields (default y direction) is a relatively complex map. To accommodate this, the values of a field in the positive (up) and negative (down) directions can be stored in separate fields.

Field3D f(0.0); // f allocated, set to zero

f.yup() // error; f.yup not allocated

f.clearParallelSlices(); // f.yup_fields and f.ydown_fields are now empty
f.yup() // error; f.yup not allocated

To have separate fields for yup and ydown, first call

f.splitParallelSlices(); // f.yup() and f.ydown() separate

f.yup(); // ok
f.yup()(0,1,0) // error; f.yup not allocated

f.yup() = 1.0; // Set f.yup() field to 1.0

f.yup()(0,1,0) // ok

Unnamed Group

Field3D &operator=(const Field3D &rhs)#

Assignment operators

Field3D &operator=(Field3D &&rhs)#
Field3D &operator=(const Field2D &rhs)#
void operator=(const FieldPerp &rhs)#

return void, as only part initialised

Field3D &operator=(BoutReal val)#

Unnamed Group

Field3D &operator+=(const Field3D &rhs)#

Addition operators

Field3D &operator+=(const Field2D &rhs)#
Field3D &operator+=(BoutReal rhs)#

Unnamed Group

Field3D &operator-=(const Field3D &rhs)#

Subtraction operators

Field3D &operator-=(const Field2D &rhs)#
Field3D &operator-=(BoutReal rhs)#

Unnamed Group

Field3D &operator*=(const Field3D &rhs)#

Multiplication operators

Field3D &operator*=(const Field2D &rhs)#
Field3D &operator*=(BoutReal rhs)#

Unnamed Group

Field3D &operator/=(const Field3D &rhs)#

Division operators

Field3D &operator/=(const Field2D &rhs)#
Field3D &operator/=(BoutReal rhs)#

Public Types

using ind_type = Ind3D#
using value_type = BoutReal#

Data type stored in this field.

Public Functions

Field3D(Mesh *localmesh = nullptr, CELL_LOC location_in = CELL_CENTRE, DirectionTypes directions_in = {YDirectionType::Standard, ZDirectionType::Standard})#

Constructor.

Constructor

Note: the global “mesh” can’t be passed here because fields may be created before the mesh is.

Field3D(const Field3D &f)#

Copy constructor

Doesn’t copy any data, just create a new reference to the same data (copy on change later)

inline Field3D(Field3D &&f) noexcept#

Move constructor.

Field3D(const Field2D &f)#

Constructor from 2D field.

Field3D(BoutReal val, Mesh *localmesh = nullptr)#

Constructor from value.

Field3D(Array<BoutReal> data, Mesh *localmesh, CELL_LOC location = CELL_CENTRE, DirectionTypes directions_in = {YDirectionType::Standard, ZDirectionType::Standard})#

Constructor from Array and Mesh.

~Field3D() override#

Destructor.

Field3D &allocate()#

Ensures that memory is allocated and unique

inline bool isAllocated() const#

Test if data is allocated

BOUT_HOST_DEVICE Field3D * timeDeriv ()

Return a pointer to the time-derivative field

The first time this is called, a new field will be allocated. Subsequent calls return the same field

inline virtual int getNx() const override#

Return the number of nx points

inline virtual int getNy() const override#

Return the number of ny points

inline virtual int getNz() const override#

Return the number of nz points

inline virtual Field3D &setLocation(CELL_LOC new_location) override#

Set variable location for staggered grids to

Parameters:

new_location – Throws BoutException if new_location is not CELL_CENTRE and staggered grids are turned off and checks are on. If checks are off, silently sets location to CELL_CENTRE instead.

inline virtual Field3D &setDirectionY(YDirectionType d) override#
void splitParallelSlices()#

Ensure that this field has separate fields for yup and ydown.

void clearParallelSlices()#

Clear the parallel slices, yup and ydown

inline bool hasParallelSlices() const#

Check if this field has yup and ydown fields.

inline Field3D &yup(std::vector<Field3D>::size_type index = 0)#

Check if this field has yup and ydown fields Return reference to yup field

inline const Field3D &yup(std::vector<Field3D>::size_type index = 0) const#

Return const reference to yup field.

inline Field3D &ydown(std::vector<Field3D>::size_type index = 0)#

Return reference to ydown field.

inline const Field3D &ydown(std::vector<Field3D>::size_type index = 0) const#

Return const reference to ydown field.

Field3D &ynext(int offset)#

Return the parallel slice at offset

offset of 0 returns the main field itself

const Field3D &ynext(int offset) const#
bool requiresTwistShift(bool twist_shift_enabled)#

If twist_shift_enabled is true, does this Field3D require a twist-shift at branch cuts on closed field lines?

const Region<Ind3D> &getRegion(REGION region) const#

Return a Region<Ind3D> reference to use to iterate over this field

Example

This loops over the interior points, not the boundary and inside the loop the index is used to calculate the difference between the point one index up in x (i.xp()) and one index down in x (i.xm()), putting the result into a different field ‘g’

for(const auto &i : f.getRegion(RGN_NOBNDRY)) { g[i] = f[i.xp()] - f[i.xm()]; }

const Region<Ind3D> &getRegion(const std::string &region_name) const#
const Region<Ind3D> &getValidRegionWithDefault(const std::string &region_name) const#

Use region provided by the default, and if none is set, use the provided one.

void setRegion(const std::string &region_name)#
inline void resetRegion()#
inline void setRegion(size_t id)#
inline void setRegion(std::optional<size_t> id)#
inline std::optional<size_t> getRegionID() const#
const Region<Ind2D> &getRegion2D(REGION region) const#

Return a Region<Ind2D> reference to use to iterate over the x- and y-indices of this field

const Region<Ind2D> &getRegion2D(const std::string &region_name) const#
inline Region<Ind3D>::RegionIndices::const_iterator begin() const#
inline Region<Ind3D>::RegionIndices::const_iterator end() const#
inline BoutReal &BOUT_HOST_DEVICE operator[] (const Ind3D &d)
inline const BoutReal &BOUT_HOST_DEVICE operator[] (const Ind3D &d) const
BoutReal &BOUT_HOST_DEVICE operator() (const IndPerp &d, int jy)
const BoutReal &BOUT_HOST_DEVICE operator() (const IndPerp &d, int jy) const
BoutReal &BOUT_HOST_DEVICE operator() (const Ind2D &d, int jz)
const BoutReal &BOUT_HOST_DEVICE operator() (const Ind2D &d, int jz) const
inline BoutReal &operator()(int jx, int jy, int jz)#

Direct access to the underlying data array

If CHECK > 2 then bounds checking is performed

If CHECK <= 2 then no checks are performed, to allow inlining and optimisation of inner loops

inline const BoutReal &operator()(int jx, int jy, int jz) const#
inline const BoutReal *operator()(int jx, int jy) const#

Direct access to the underlying data array

This version returns a pointer to a data array, and is intended for use with FFT routines. The data is guaranteed to be contiguous in Z index

inline BoutReal *operator()(int jx, int jy)#
inline virtual bool is3D() const override#

True if variable is 3D.

inline virtual void doneComms() override#
Field3D &calcParallelSlices()#
virtual void applyBoundary(bool init = false) override#
void applyBoundary(BoutReal t)#
void applyBoundary(const std::string &condition)#
inline void applyBoundary(const char *condition)#
void applyBoundary(const std::string &region, const std::string &condition)#
virtual void applyTDerivBoundary() override#
void setBoundaryTo(const Field3D &f3d)#

Copy the boundary values half-way between cells This uses 2nd order central differences to set the value on the boundary to the value on the boundary in field f3d. Note: does not just copy values in boundary region.

virtual void applyParallelBoundary() override#
virtual void applyParallelBoundary(BoutReal t) override#
virtual void applyParallelBoundary(const std::string &condition) override#
virtual void applyParallelBoundary(const std::string &region, const std::string &condition) override#
void applyParallelBoundary(const std::string &region, const std::string &condition, Field3D *f)#
inline virtual int size() const override#

Get the total number of points.

Private Members

int nx = {-1}#

Array sizes (from fieldmesh). These are valid only if fieldmesh is not null.

int ny = {-1}#
int nz = {-1}#
Array<BoutReal> data#

Internal data array. Handles allocation/freeing of memory.

Field3D *deriv = {nullptr}#

Time derivative (may be nullptr)

std::vector<Field3D> yup_fields = {}#

Fields containing values along Y.

std::vector<Field3D> ydown_fields = {}#
std::optional<size_t> regionID#

RegionID over which the field is valid.

Friends

friend class Vector3D
friend class Vector2D
friend void swap(Field3D &first, Field3D &second) noexcept#