Holds field data, performs FFTs and HDF5 I/O operations. More...

#include <field.hpp>

Inheritance diagram for field< rnumber, be, fc >:

Collaboration diagram for field< rnumber, be, fc >:

Public Member Functions
	field (const int nx, const int ny, const int nz, const MPI_Comm COMM_TO_USE, const unsigned FFTW_PLAN_RIGOR=DEFAULT_FFTW_FLAG)

	~field () noexcept(false)

int	io (const std::string fname, const std::string field_name, const int iteration, const bool read=true)

int	io_database (const std::string fname, const std::string field_name, const int toffset, const bool read=true)

int	write_0slice (const hid_t group, const std::string field_name, const int iteration)

int	write_filtered (const std::string fname, const std::string field_name, const int iteration, const int nx, const int ny, const int nz)

int	io_binary (const std::string fname, const int iteration, const bool read=true)

void	dft ()

void	ift ()

void	normalize ()

void	symmetrize ()
	Enforce Hermitian symmetry (fast and reasonable)

void	symmetrize_FFT ()
	Enforce Hermitian symmetry (slow, reference)

void	symmetrize_alternate ()
	Enforce Hermitian symmetry (unoptimized, amplitude-aware)

void	Hermitian_reflect ()

void	compute_rspace_xincrement_stats (const int xcells, const hid_t group, const std::string dset_name, const hsize_t toffset, const std::vector< double > max_estimate, field< rnumber, be, fc > *tmp_field=NULL)

void	compute_rspace_stats (const hid_t group, const std::string dset_name, const hsize_t toffset, const std::vector< double > max_estimate, const unsigned int maximum_moment=9)
	Compute statistics (real-space representation)

double	compute_CFL_velocity ()
	Compute CFL velocity.

void	compute_rspace_zaverage (const hid_t group, const std::string dset_name, const hsize_t toffset)

int	get_nx () const

int	get_ny () const

int	get_nz () const

rnumber *__restrict__	get_rdata ()

const rnumber *__restrict__	get_rdata () const

fftw_interface< rnumber >::complex *__restrict__	get_cdata ()

fftw_interface< rnumber >::complex *__restrict__	get_cdata () const

rnumber &	rval (ptrdiff_t rindex, unsigned int component=0)

const rnumber &	rval (ptrdiff_t rindex, unsigned int component=0) const

rnumber &	rval (ptrdiff_t rindex, int comp1, int comp0)

rnumber &	cval (ptrdiff_t cindex, int imag)

rnumber &	cval (ptrdiff_t cindex, int component, int imag)

rnumber &	cval (ptrdiff_t cindex, int comp1, int comp0, int imag)

field< rnumber, be, fc > &	operator= (const typename fftw_interface< rnumber >::complex *__restrict__ source)

field< rnumber, be, fc > &	operator= (const rnumber *__restrict__ source)

field< rnumber, be, fc > &	operator= (const rnumber value)

field< rnumber, be, fc > &	operator= (const field< rnumber, be, fc > &src)

template<kspace_dealias_type dt>
void	compute_stats (kspace< be, dt > *kk, const hid_t group, const std::string dset_name, const hsize_t toffset, const double max_estimate)
	Compute field statistics.

template<kspace_dealias_type dt>
double	L2norm (kspace< be, dt > *kk)

void	impose_zero_mode ()

template<class func_type >
void	RLOOP_simd (func_type expression)

template<class func_type >
void	RLOOP (func_type expression)

ptrdiff_t	get_cindex (ptrdiff_t xindex, ptrdiff_t yindex, ptrdiff_t zindex)

ptrdiff_t	get_rindex (ptrdiff_t xindex, ptrdiff_t yindex, ptrdiff_t zindex) const

ptrdiff_t	get_rindex_from_global (const ptrdiff_t in_global_x, const ptrdiff_t in_global_y, const ptrdiff_t in_global_z) const

int	print_plan (const std::string preamble)

Static Public Member Functions
static constexpr int	get_number_of_components (void)

Data Fields
hsize_t	npoints

bool	real_space_representation

int	myrank

int	nprocs

MPI_Comm	comm

field_layout< fc > *	clayout

field_layout< fc > *	rlayout

field_layout< fc > *	rmemlayout

fftw_interface< rnumber >::many_plan	c2r_plan

fftw_interface< rnumber >::many_plan	r2c_plan

unsigned	fftw_plan_rigor

hid_t	rnumber_H5T

hid_t	cnumber_H5T

Static Public Attributes
static constexpr int	number_of_components = ncomp(fc)

Private Attributes
rnumber *__restrict__	data

Detailed Description

template<typename rnumber, field_backend be, field_components fc>
class field< rnumber, be, fc >

Holds field data, performs FFTs and HDF5 I/O operations.

The purpose of this class is to manage memory for field data, create/destroy FFT plans for them, and compute HDF5 input/output operations.

FFTW recommendations are to create different plans for different arrays, hence the plans are member variables. All plans are for in-place transforms, since even with out-of-place transforms there are no guarantees that input data is not messed up by an inverse FFT, so there's no point in wasting the memory.

Constructor & Destructor Documentation

◆ field()

template<typename rnumber , field_backend be, field_components fc>

field< rnumber, be, fc >::field	(	const int	nx,
		const int	ny,
		const int	nz,
		const MPI_Comm	COMM_TO_USE,
		const unsigned	FFTW_PLAN_RIGOR = DEFAULT_FFTW_FLAG )

◆ ~field()

template<typename rnumber , field_backend be, field_components fc>

field< rnumber, be, fc >::~field ( )

Member Function Documentation

◆ compute_CFL_velocity()

template<typename rnumber , field_backend be, field_components fc>

double field< rnumber, be, fc >::compute_CFL_velocity ( )

Compute CFL velocity.

For the CFL condition we need a specific "maximum" norm of the velocity field.

\[ \| u \|_{CFL} \equiv \max_{\mathbf{x}}\left(\sum_{i=1,2,3}\big(|u_i(\mathbf{x})|\big)\right) \]

◆ compute_rspace_stats()

template<typename rnumber , field_backend be, field_components fc>

void field< rnumber, be, fc >::compute_rspace_stats	(	const hid_t	group,
		const std::string	dset_name,
		const hsize_t	toffset,
		const std::vector< double >	max_estimate,
		const unsigned int	maximum_moment = 9 )

Compute statistics (real-space representation)

Note

Moments from 1 through 8 are computed. Moment array will have 10 elements:

index 0 contains minimum value
indices 1 through 8 contain moments 1 through 8
index 9 contains maximum value

◆ compute_rspace_xincrement_stats()

template<typename rnumber , field_backend be, field_components fc>

void field< rnumber, be, fc >::compute_rspace_xincrement_stats	(	const int	xcells,
		const hid_t	group,
		const std::string	dset_name,
		const hsize_t	toffset,
		const std::vector< double >	max_estimate,
		field< rnumber, be, fc > *	tmp_field = NULL )

◆ compute_rspace_zaverage()

template<typename rnumber , field_backend be, field_components fc>

void field< rnumber, be, fc >::compute_rspace_zaverage	(	const hid_t	group,
		const std::string	dset_name,
		const hsize_t	toffset )

◆ compute_stats()

template<typename rnumber , field_backend be, field_components fc>

template<kspace_dealias_type dt>

template void field< rnumber, be, fc >::compute_stats< SMOOTH >	(	kspace< be, dt > *	kk,
		const hid_t	group,
		const std::string	dset_name,
		const hsize_t	toffset,
		const double	max_estimate )

Compute field statistics.

Warning: This method calls compute_rspace_stats, as well as cospectrum This means that there will be a Fourier transform performed at some point.; The "max_estimate" value is manipulated before the call to compute_rspace_stats: if the field is a vector field, the maximum estimate for the field magnitude is multiplied by sqrt(3). This is a naive attempt to account for the fact that the magnitude is larger than any 1 of the 3 components.

◆ cval() [1/3]

template<typename rnumber , field_backend be, field_components fc>

rnumber & field< rnumber, be, fc >::cval	(	ptrdiff_t	cindex,
		int	comp1,
		int	comp0,
		int	imag )

inline

◆ cval() [2/3]

template<typename rnumber , field_backend be, field_components fc>

rnumber & field< rnumber, be, fc >::cval	(	ptrdiff_t	cindex,
		int	component,
		int	imag )

inline

◆ cval() [3/3]

template<typename rnumber , field_backend be, field_components fc>

rnumber & field< rnumber, be, fc >::cval	(	ptrdiff_t	cindex,
		int	imag )

inline

◆ dft()

template<typename rnumber , field_backend be, field_components fc>

void field< rnumber, be, fc >::dft ( )

◆ get_cdata() [1/2]

template<typename rnumber , field_backend be, field_components fc>

fftw_interface< rnumber >::complex *__restrict__ field< rnumber, be, fc >::get_cdata ( )

inline

◆ get_cdata() [2/2]

template<typename rnumber , field_backend be, field_components fc>

fftw_interface< rnumber >::complex *__restrict__ field< rnumber, be, fc >::get_cdata ( ) const

inline

◆ get_cindex()

template<typename rnumber , field_backend be, field_components fc>

ptrdiff_t field< rnumber, be, fc >::get_cindex	(	ptrdiff_t	xindex,
		ptrdiff_t	yindex,
		ptrdiff_t	zindex )

inline

◆ get_number_of_components()

template<typename rnumber , field_backend be, field_components fc>

static constexpr int field< rnumber, be, fc >::get_number_of_components ( void )

inlinestaticconstexpr

◆ get_nx()

template<typename rnumber , field_backend be, field_components fc>

int field< rnumber, be, fc >::get_nx ( ) const

inline

◆ get_ny()

template<typename rnumber , field_backend be, field_components fc>

int field< rnumber, be, fc >::get_ny ( ) const

inline

◆ get_nz()

template<typename rnumber , field_backend be, field_components fc>

int field< rnumber, be, fc >::get_nz ( ) const

inline

◆ get_rdata() [1/2]

template<typename rnumber , field_backend be, field_components fc>

rnumber *__restrict__ field< rnumber, be, fc >::get_rdata ( )

inline

◆ get_rdata() [2/2]

template<typename rnumber , field_backend be, field_components fc>

const rnumber *__restrict__ field< rnumber, be, fc >::get_rdata ( ) const

inline

◆ get_rindex()

template<typename rnumber , field_backend be, field_components fc>

ptrdiff_t field< rnumber, be, fc >::get_rindex	(	ptrdiff_t	xindex,
		ptrdiff_t	yindex,
		ptrdiff_t	zindex ) const

inline

◆ get_rindex_from_global()

template<typename rnumber , field_backend be, field_components fc>

ptrdiff_t field< rnumber, be, fc >::get_rindex_from_global	(	const ptrdiff_t	in_global_x,
		const ptrdiff_t	in_global_y,
		const ptrdiff_t	in_global_z ) const

inline

◆ Hermitian_reflect()

template<typename rnumber , field_backend be, field_components fc>

void field< rnumber, be, fc >::Hermitian_reflect ( )

◆ ift()

template<typename rnumber , field_backend be, field_components fc>

void field< rnumber, be, fc >::ift ( )

◆ impose_zero_mode()

template<typename rnumber , field_backend be, field_components fc>

void field< rnumber, be, fc >::impose_zero_mode ( )

inline

◆ io()

template<typename rnumber , field_backend be, field_components fc>

int field< rnumber, be, fc >::io	(	const std::string	fname,
		const std::string	field_name,
		const int	iteration,
		const bool	read = true )

◆ io_binary()

template<typename rnumber , field_backend be, field_components fc>

int field< rnumber, be, fc >::io_binary	(	const std::string	fname,
		const int	iteration,
		const bool	read = true )

◆ io_database()

template<typename rnumber , field_backend be, field_components fc>

int field< rnumber, be, fc >::io_database	(	const std::string	fname,
		const std::string	field_name,
		const int	toffset,
		const bool	read = true )

◆ L2norm()

template<typename rnumber , field_backend be, field_components fc>

template<kspace_dealias_type dt>

template double field< rnumber, be, fc >::L2norm< SMOOTH > ( kspace< be, dt > * kk )

◆ normalize()

template<typename rnumber , field_backend be, field_components fc>

void field< rnumber, be, fc >::normalize ( )

◆ operator=() [1/4]

template<typename rnumber , field_backend be, field_components fc>

field< rnumber, be, fc > & field< rnumber, be, fc >::operator= ( const field< rnumber, be, fc > & src )

◆ operator=() [2/4]

template<typename rnumber , field_backend be, field_components fc>

field< rnumber, be, fc > & field< rnumber, be, fc >::operator= ( const rnumber *__restrict__ source )

◆ operator=() [3/4]

template<typename rnumber , field_backend be, field_components fc>

field< rnumber, be, fc > & field< rnumber, be, fc >::operator= ( const rnumber value )

◆ operator=() [4/4]

template<typename rnumber , field_backend be, field_components fc>

field< rnumber, be, fc > & field< rnumber, be, fc >::operator= ( const typename fftw_interface< rnumber >::complex *__restrict__ source )

◆ print_plan()

template<typename rnumber , field_backend be, field_components fc>

int field< rnumber, be, fc >::print_plan ( const std::string preamble )

◆ RLOOP()

template<typename rnumber , field_backend be, field_components fc>

template<class func_type >

void field< rnumber, be, fc >::RLOOP ( func_type expression )

inline

◆ RLOOP_simd()

template<typename rnumber , field_backend be, field_components fc>

template<class func_type >

void field< rnumber, be, fc >::RLOOP_simd ( func_type expression )

inline

◆ rval() [1/3]

template<typename rnumber , field_backend be, field_components fc>

rnumber & field< rnumber, be, fc >::rval	(	ptrdiff_t	rindex,
		int	comp1,
		int	comp0 )

inline

◆ rval() [2/3]

template<typename rnumber , field_backend be, field_components fc>

rnumber & field< rnumber, be, fc >::rval	(	ptrdiff_t	rindex,
		unsigned int	component = 0 )

inline

◆ rval() [3/3]

template<typename rnumber , field_backend be, field_components fc>

const rnumber & field< rnumber, be, fc >::rval	(	ptrdiff_t	rindex,
		unsigned int	component = 0 ) const

inline

◆ symmetrize()

template<typename rnumber , field_backend be, field_components fc>

void field< rnumber, be, fc >::symmetrize ( )

Enforce Hermitian symmetry (fast and reasonable)

TurTLE uses real-to-complex and complex-to-real FFTW transforms, because the equations are PDEs of real-valued fields. Hermitian symmetry means that Fourier-transformed real valued fields must respect the equation \(\hat f(-\mathbf{k}) = {\hat f}^*(\mathbf{k})\).

FFTW enforces this property mainly by only storing the positive half of the Fourier grid for the fastest array component. In TurTLE's case, this means \( k_x \). For the \( k_x = 0 \) plane, the symmetry must be enforced.

This method uses an arithmetic mean of the \( (0, k_y, k_z)\) mode and the conjugate of the \((0, -k_y, -k_z) \) mode to generate the desired values. The method is fast (other than the required MPI communications).

Note: the method is adequate either in cases where deviations from Hermitian symmetry are small, or in cases where deviations from correct physics is irrelevant. In practice: initial condition fields may be strongly perturbed by the application of this method, but they are unphysical anyway; during the quasistationary regime of some simulation, the method is applied regularly to all relevant fields, and deviations are expected to be small, i.e. effect on PDE approximations is negligible.

◆ symmetrize_alternate()

template<typename rnumber , field_backend be, field_components fc>

void field< rnumber, be, fc >::symmetrize_alternate ( )

Enforce Hermitian symmetry (unoptimized, amplitude-aware)

TurTLE uses real-to-complex and complex-to-real FFTW transforms, because the equations are PDEs of real-valued fields. Hermitian symmetry means that Fourier-transformed real valued fields must respect the equation \(\hat f(-\mathbf{k}) = {\hat f}^*(\mathbf{k})\).

FFTW enforces this property mainly by only storing the positive half of the Fourier grid for the fastest array component. In TurTLE's case, this means \( k_x \). For the \( k_x = 0 \) plane, the symmetry must be enforced.

This method is an alternative to the default arithmetic mean method, meant to be used in special circumstances where it is important to retain the exact amplitude of modes. Rather than an arithmetic mean, here we compute the amplitudes and phases for the \((0, k_y, k_z)\) and \((0, -k_y, -k_z)\) modes. We then compute a mean amplitude as the square root of the product of the two amplitudes, and a mean phase as the arithmetic mean of the two phases. When this method is applied to a field with fixed amplitudes, but random phases, it should preserve the spectrum of the initial field exactly.

◆ symmetrize_FFT()

template<typename rnumber , field_backend be, field_components fc>

void field< rnumber, be, fc >::symmetrize_FFT ( )

Enforce Hermitian symmetry (slow, reference)

TurTLE uses real-to-complex and complex-to-real FFTW transforms, because the equations are PDEs of real-valued fields. Hermitian symmetry means that Fourier-transformed real valued fields must respect the equation \(\hat f(-\mathbf{k}) = {\hat f}^*(\mathbf{k})\).

FFTW enforces this property mainly by only storing the positive half of the Fourier grid for the fastest array component. In TurTLE's case, this means \( k_x \). For the \( k_x = 0 \) plane, the symmetry must be enforced.

This method uses a pair of backwards-forwards FFTs, which leads to FFTW effectively imposing Hermitian symmetry. It should be used as a reference implementation to calibrate against in the general case.

◆ write_0slice()

template<typename rnumber , field_backend be, field_components fc>

int field< rnumber, be, fc >::write_0slice	(	const hid_t	group,
		const std::string	field_name,
		const int	iteration )

◆ write_filtered()

template<typename rnumber , field_backend be, field_components fc>

int field< rnumber, be, fc >::write_filtered	(	const std::string	fname,
		const std::string	field_name,
		const int	iteration,
		const int	nx,
		const int	ny,
		const int	nz )

set up counts and offsets x is easy, since only positive modes are present

three options for y: this->starts[0] <= ny/2 ny / 2 < this->starts[0] +this->clayout->subsizes[0] < this->sizes[0] - ny/2 this->starts[0] >= this->sizes[0] - ny/2 we don't care about saving the ny/2 mode, because of symmetry

for z, we need to take into account that there are both positive and negative modes

Field Documentation

◆ c2r_plan

template<typename rnumber , field_backend be, field_components fc>

fftw_interface<rnumber>::many_plan field< rnumber, be, fc >::c2r_plan

◆ clayout

template<typename rnumber , field_backend be, field_components fc>

field_layout<fc>* field< rnumber, be, fc >::clayout

◆ cnumber_H5T

template<typename rnumber , field_backend be, field_components fc>

hid_t field< rnumber, be, fc >::cnumber_H5T

◆ comm

template<typename rnumber , field_backend be, field_components fc>

MPI_Comm field< rnumber, be, fc >::comm

MPI communicator this fields lives in.

◆ data

template<typename rnumber , field_backend be, field_components fc>

rnumber* __restrict__ field< rnumber, be, fc >::data

private

data array

◆ fftw_plan_rigor

template<typename rnumber , field_backend be, field_components fc>

unsigned field< rnumber, be, fc >::fftw_plan_rigor

◆ myrank

template<typename rnumber , field_backend be, field_components fc>

int field< rnumber, be, fc >::myrank

◆ npoints

template<typename rnumber , field_backend be, field_components fc>

hsize_t field< rnumber, be, fc >::npoints

total number of grid points. Useful for normalization.

◆ nprocs

template<typename rnumber , field_backend be, field_components fc>

int field< rnumber, be, fc >::nprocs

basic MPI information.

◆ number_of_components

template<typename rnumber , field_backend be, field_components fc>

constexpr int field< rnumber, be, fc >::number_of_components = ncomp(fc)

staticconstexpr

◆ r2c_plan

template<typename rnumber , field_backend be, field_components fc>

fftw_interface<rnumber>::many_plan field< rnumber, be, fc >::r2c_plan

◆ real_space_representation

template<typename rnumber , field_backend be, field_components fc>

bool field< rnumber, be, fc >::real_space_representation

true if field is in real space representation.

◆ rlayout

template<typename rnumber , field_backend be, field_components fc>

field_layout<fc> * field< rnumber, be, fc >::rlayout

◆ rmemlayout

template<typename rnumber , field_backend be, field_components fc>

field_layout<fc> * field< rnumber, be, fc >::rmemlayout

◆ rnumber_H5T

template<typename rnumber , field_backend be, field_components fc>

hid_t field< rnumber, be, fc >::rnumber_H5T

The documentation for this class was generated from the following files:

/builds/TurTLE/turtle/cpp/field.hpp
/builds/TurTLE/turtle/cpp/field.cpp

Public Member Functions

Static Public Member Functions

Data Fields

Static Public Attributes

Private Attributes

Detailed Description

Constructor & Destructor Documentation

◆ field()

◆ ~field()

Member Function Documentation

◆ compute_CFL_velocity()

◆ compute_rspace_stats()

◆ compute_rspace_xincrement_stats()

◆ compute_rspace_zaverage()

◆ compute_stats()

◆ cval() [1/3]

◆ cval() [2/3]

◆ cval() [3/3]

◆ dft()

◆ get_cdata() [1/2]

◆ get_cdata() [2/2]

◆ get_cindex()

◆ get_number_of_components()

◆ get_nx()

◆ get_ny()

◆ get_nz()

◆ get_rdata() [1/2]

◆ get_rdata() [2/2]

◆ get_rindex()

◆ get_rindex_from_global()

◆ Hermitian_reflect()

◆ ift()

◆ impose_zero_mode()

◆ io()

◆ io_binary()

◆ io_database()

◆ L2norm()

◆ normalize()

◆ operator=() [1/4]

◆ operator=() [2/4]

◆ operator=() [3/4]

◆ operator=() [4/4]

◆ print_plan()

◆ RLOOP()

◆ RLOOP_simd()

◆ rval() [1/3]

◆ rval() [2/3]

◆ rval() [3/3]

◆ symmetrize()

◆ symmetrize_alternate()

◆ symmetrize_FFT()

◆ write_0slice()

◆ write_filtered()

Field Documentation

◆ c2r_plan

◆ clayout

◆ cnumber_H5T

◆ comm

◆ data

◆ fftw_plan_rigor

◆ myrank

◆ npoints

◆ nprocs

◆ number_of_components

◆ r2c_plan

◆ real_space_representation

◆ rlayout

◆ rmemlayout

◆ rnumber_H5T