JDFTx  1.7.0
Operator implementations

Files

file  BlasExtra.h
 Commonly used BLAS-like routines.
 
file  Operators.h
 Operators on ScalarField's and ScalarFieldTilde's.
 
file  TranslationOperator.h
 

Classes

class  ColumnBundleTransform
 Handle transformation of ColumnBundles upon symmetry operations. More...
 
class  TranslationOperator
 Abstract base class for translation operators. More...
 
class  TranslationOperatorSpline
 Translation operator which works in real space using interpolating splines. More...
 
class  TranslationOperatorFourier
 The exact translation operator in PW basis, although much slower and with potential ringing issues. More...
 

Macros

#define callPref(functionName)   functionName##_gpu
 Select between functionName and functionName_gpu for the CPU and GPU executables respectively.
 
#define Tptr   std::shared_ptr<T>
 shorthand for writing the template operators (undef'd at end of header)
 

Functions

template<typename Ty , typename Tx >
void eblas_mul (const int N, const Tx *X, const int incX, Ty *Y, const int incY)
 Templated elementwise multiply Y *= X for arrays X, Y. More...
 
void eblas_dmul (const int N, const double *X, const int incX, double *Y, const int incY)
 Specialization of eblas_mul() for double[] *= double[].
 
void eblas_zmul (const int N, const complex *X, const int incX, complex *Y, const int incY)
 Specialization of eblas_mul() for complex[] *= complex[].
 
void eblas_zmuld (const int N, const double *X, const int incX, complex *Y, const int incY)
 Specialization of eblas_mul() for complex[] *= double[].
 
void eblas_dmul_gpu (const int N, const double *X, const int incX, double *Y, const int incY)
 Equivalent of eblas_dmul() for GPU data pointers.
 
void eblas_zmul_gpu (const int N, const complex *X, const int incX, complex *Y, const int incY)
 Equivalent of eblas_zmul() for GPU data pointers.
 
void eblas_zmuld_gpu (const int N, const double *X, const int incX, complex *Y, const int incY)
 Equivalent of eblas_zmuld() for GPU data pointers.
 
template<typename Ty , typename Tx >
void eblas_div (const int N, const Tx *X, const int incX, Ty *Y, const int incY)
 Templated elementwise divide Y /= X for arrays X, Y. More...
 
void eblas_ddiv (const int N, const double *X, const int incX, double *Y, const int incY)
 Specialization of eblas_div() for double[] /= double[].
 
void eblas_zdiv (const int N, const complex *X, const int incX, complex *Y, const int incY)
 Specialization of eblas_div() for complex[] /= complex[].
 
void eblas_zdivd (const int N, const double *X, const int incX, complex *Y, const int incY)
 Specialization of eblas_div() for complex[] /= double[].
 
void eblas_ddiv_gpu (const int N, const double *X, const int incX, double *Y, const int incY)
 Equivalent of eblas_ddiv() for GPU data pointers.
 
void eblas_zdiv_gpu (const int N, const complex *X, const int incX, complex *Y, const int incY)
 Equivalent of eblas_zdiv() for GPU data pointers.
 
void eblas_zdivd_gpu (const int N, const double *X, const int incX, complex *Y, const int incY)
 Equivalent of eblas_zdivd() for GPU data pointers.
 
void eblas_sumStrided (const int N, const int stride, const double *X, double *result)
 Return sum over array of doubles with specified stride.
 
void eblas_sumStrided_gpu (const int N, const int stride, const double *X, double *result)
 Equivalent of eblas_sumStrided() for X on GPU; result is on CPU.
 
template<typename T >
eblas_sum (const int N, const T *X)
 Return sum over array of generic type T composed only of a certain number of double members.
 
template<typename T >
eblas_sum_gpu (const int N, const T *X)
 Equivalent of eblas_sum() for GPU data pointers.
 
void eblas_lincomb (const int N, const complex &sX, const complex *X, const int incX, const complex &sY, const complex *Y, const int incY, complex *Z, const int incZ)
 Elementwise linear combination Z = sX * X + sY * Y. More...
 
void eblas_lincomb_gpu (const int N, const complex &sX, const complex *X, const int incX, const complex &sY, const complex *Y, const int incY, complex *Z, const int incZ)
 Equivalent of eblas_lincomb() for GPU data pointers.
 
void eblas_zgemm (CBLAS_TRANSPOSE TransA, CBLAS_TRANSPOSE TransB, int M, int N, int K, const complex &alpha, const complex *A, const int lda, const complex *B, const int ldb, const complex &beta, complex *C, const int ldc)
 Threaded complex matrix multiply (threaded wrapper around zgemm) All the parameters have the same meaning as in cblas_zgemm, except element order is always Column Major (FORTRAN order!)
 
void eblas_zgemm_gpu (CBLAS_TRANSPOSE TransA, CBLAS_TRANSPOSE TransB, int M, int N, int K, const complex &alpha, const complex *A, const int lda, const complex *B, const int ldb, const complex &beta, complex *C, const int ldc)
 Wrap cublasZgemm to provide the same interface as eblas_zgemm()
 
void eblas_scatter_zdaxpy (const int Nindex, double a, const int *index, const complex *x, complex *y, bool conjx=false, const complex *w=0, bool conjw=false)
 Scatter y(index) += a * x. More...
 
void eblas_scatter_zaxpy (const int Nindex, complex a, const int *index, const complex *x, complex *y, bool conjx=false, const complex *w=0, bool conjw=false)
 Equivalent of eblas_scatter_zdaxpy() with a complex scale factor.
 
void eblas_scatter_daxpy (const int Nindex, double a, const int *index, const double *x, double *y, const double *w=0)
 Equivalent of eblas_scatter_zdaxpy() for real data arrays.
 
void eblas_gather_zdaxpy (const int Nindex, double a, const int *index, const complex *x, complex *y, bool conjx=false, const complex *w=0, bool conjw=false)
 Gather y += a * x(index) More...
 
void eblas_gather_zaxpy (const int Nindex, complex a, const int *index, const complex *x, complex *y, bool conjx=false, const complex *w=0, bool conjw=false)
 Equivalent of eblas_gather_zdaxpy() with a complex scale factor.
 
void eblas_gather_daxpy (const int Nindex, double a, const int *index, const double *x, double *y, const double *w=0)
 Equivalent of eblas_scatter_zdaxpy() for real data arrays.
 
void eblas_scatter_zdaxpy_gpu (const int Nindex, double a, const int *index, const complex *x, complex *y, bool conjx=false, const complex *w=0, bool conjw=false)
 Equivalent of eblas_scatter_zdaxpy() for GPU data pointers.
 
void eblas_scatter_zaxpy_gpu (const int Nindex, complex a, const int *index, const complex *x, complex *y, bool conjx=false, const complex *w=0, bool conjw=false)
 Equivalent of eblas_scatter_zaxpy() for GPU data pointers.
 
void eblas_scatter_daxpy_gpu (const int Nindex, double a, const int *index, const double *x, double *y, const double *w=0)
 Equivalent of eblas_scatter_daxpy() for GPU data pointers.
 
void eblas_gather_zdaxpy_gpu (const int Nindex, double a, const int *index, const complex *x, complex *y, bool conjx=false, const complex *w=0, bool conjw=false)
 Equivalent of eblas_gather_zdaxpy() for GPU data pointers.
 
void eblas_gather_zaxpy_gpu (const int Nindex, complex a, const int *index, const complex *x, complex *y, bool conjx=false, const complex *w=0, bool conjw=false)
 Equivalent of eblas_gather_zaxpy() for GPU data pointers.
 
void eblas_gather_daxpy_gpu (const int Nindex, double a, const int *index, const double *x, double *y, const double *w=0)
 Equivalent of eblas_gather_daxpy() for GPU data pointers.
 
void eblas_accumNorm (int N, const double &a, const complex *x, double *y)
 Accumulate elementwise norm of a complex array x into y i.e. y += a x conj(x) More...
 
void eblas_accumProd (int N, const double &a, const complex *xU, const complex *xC, double *yRe, double *yIm)
 Accumulate elementwise product of two complex arrays xU and xC into real and imaginary parts yRe and yIm i.e. (yRe + i yIm) += a xU conj(xC) More...
 
void eblas_accumProdComplex (int N, const double &a, const complex *xU, const complex *xC, complex *y)
 Accumulate elementwise product of two complex arrays xU and xC into y. More...
 
void eblas_accumNorm_gpu (int N, const double &a, const complex *x, double *y)
 Equivalent of eblas_accumNorm() for GPU data pointers.
 
void eblas_accumProd_gpu (int N, const double &a, const complex *xU, const complex *xC, double *yRe, double *yIm)
 Equivalent of eblas_accumProd() for GPU data pointers.
 
void eblas_accumProdComplex_gpu (int N, const double &a, const complex *xU, const complex *xC, complex *y)
 Equivalent of eblas_accumProd() for GPU data pointers.
 
void eblas_symmetrize (int N, int n, const int *symmIndex, double *x)
 Symmetrize an array x, using N n-fold equivalence classes in symmIndex. More...
 
void eblas_symmetrize (int N, int n, const int *symmIndex, complex *x)
 Equivalent of eblas_symmetrize() for complex data pointers.
 
void eblas_symmetrize_gpu (int N, int n, const int *symmIndex, double *x)
 Equivalent of eblas_symmetrize() for real GPU data pointers.
 
void eblas_symmetrize_gpu (int N, int n, const int *symmIndex, complex *x)
 Equivalent of eblas_symmetrize() for complex GPU data pointers.
 
void eblas_symmetrize (int N, int n, const int *symmIndex, const int *symmMult, const complex *phase, complex *x)
 Symmetrize a complex array x with phase factors, using N n-fold equivalence classes in symmIndex (useful for space group symmetrization in reciprocal space) More...
 
void eblas_symmetrize_gpu (int N, int n, const int *symmIndex, const int *symmMult, const complex *phase, complex *x)
 Equivalent of eblas_symmetrize() for complex GPU data pointers.
 
void eblas_symmetrize (int N, int n, const int *symmIndex, const int *symmMult, const complex *phase, const matrix3<> *rotSpin, std::vector< complex * > x)
 Symmetrize a quadruplet of complex arrays with phase factors, using N n-fold equivalence classes in symmIndex (useful for space group symmetrization of spin density matrices in reciprocal space) More...
 
void eblas_symmetrize_gpu (int N, int n, const int *symmIndex, const int *symmMult, const complex *phase, const matrix3<> *rotSpin, std::vector< complex * > x)
 Equivalent of eblas_symmetrize() for complex GPU data pointers.
 
template<typename T >
void eblas_copy (T *dest, const T *src, int N)
 Copy a data array. More...
 
template<typename T >
void eblas_zero (int N, T *x)
 Zero a data array. More...
 
void eblas_dscal (int N, double a, double *x, int incx)
 Scale a real array: threaded wrapper to the cblas_dscal BLAS1 function.
 
void eblas_zdscal (int N, double a, complex *x, int incx)
 Scale a complex array by a real scale factor: threaded wrapper to the cblas_zdscal BLAS1 function.
 
void eblas_zscal (int N, const complex &a, complex *x, int incx)
 Scale a complex array by a complex scale factor: threaded wrapper to the cblas_zscal BLAS1 function.
 
void eblas_daxpy (int N, double a, const double *x, int incx, double *y, int incy)
 Scaled-accumulate on real arrays: threaded wrapper to the cblas_daxpy BLAS1 function.
 
void eblas_zaxpy (int N, const complex &a, const complex *x, int incx, complex *y, int incy)
 Scaled-accumulate on complex arrays: threaded wrapper to the cblas_zaxpy BLAS1 function.
 
complex eblas_zdotc (int N, const complex *x, int incx, const complex *y, int incy)
 Dot product of complex arrays: threaded wrapper to the cblas_zdotc BLAS1 function.
 
double eblas_ddot (int N, const double *x, int incx, const double *y, int ncy)
 Dot product of real arrays: threaded wrapper to the cblas_ddot BLAS1 function.
 
double eblas_dznrm2 (int N, const complex *x, int incx)
 2-norm of a complex array: threaded wrapper to the cblas_dznrm2 BLAS1 function
 
double eblas_dnrm2 (int N, const double *x, int incx)
 2-norm of a real array: threaded wrapper to the cblas_dnrm2 BLAS1 function
 
template<typename T >
void eblas_copy_gpu (T *dest, const T *src, int N)
 Equivalent of eblas_copy() for GPU data pointers.
 
template<typename T >
void eblas_zero_gpu (int N, T *x)
 Equivalent of eblas_zero() for GPU data pointers.
 
void eblas_dscal_gpu (int N, double a, double *x, int incx)
 Equivalent of eblas_dscal() for GPU data pointers.
 
void eblas_zdscal_gpu (int N, double a, complex *x, int incx)
 Equivalent of eblas_zdscal() for GPU data pointers.
 
void eblas_zscal_gpu (int N, const complex &a, complex *x, int incx)
 Equivalent of eblas_zscal for GPU data pointers.
 
void eblas_daxpy_gpu (int N, double a, const double *x, int incx, double *y, int incy)
 Equivalent of eblas_daxpy() for GPU data pointers.
 
void eblas_zaxpy_gpu (int N, const complex &a, const complex *x, int incx, complex *y, int incy)
 Equivalent of eblas_zaxpy() for GPU data pointers.
 
complex eblas_zdotc_gpu (int N, const complex *x, int incx, const complex *y, int incy)
 Equivalent of eblas_zdotc() for GPU data pointers.
 
double eblas_ddot_gpu (int N, const double *x, int incx, const double *y, int incy)
 Equivalent of eblas_ddot() for GPU data pointers.
 
double eblas_dznrm2_gpu (int N, const complex *x, int incx)
 Equivalent of eblas_dznrm2() for GPU data pointers.
 
double eblas_dnrm2_gpu (int N, const double *x, int incx)
 Equivalent of eblas_dnrm2() for GPU data pointers.
 
void eblas_capMinMax (const int N, double *x, double &xMin, double &xMax, double capLo=-DBL_MAX, double capHi=+DBL_MAX)
 Find the minimum and maximum of a data array and optionally cap it from above and/or below. More...
 
void eblas_capMinMax_gpu (const int N, double *x, double &xMin, double &xMax, double capLo=-DBL_MAX, double capHi=+DBL_MAX)
 Equivalent of eblas_capMinMax() for GPU data pointers.
 
ScalarField Real (const complexScalarField &)
 real part of a complex scalar field (real-space)
 
ScalarFieldTilde Real (const complexScalarFieldTilde &)
 real part of a complex scalar field (reciprocal space)
 
ScalarField Imag (const complexScalarField &)
 imaginary part of a complex scalar field (real-space)
 
ScalarFieldTilde Imag (const complexScalarFieldTilde &)
 imaginary part of a complex scalar field (reciprocal space)
 
complexScalarField Complex (const ScalarField &)
 convert real to complex scalar field with zero imaginary part (real-space)
 
complexScalarField Complex (const ScalarField &re, const ScalarField &im)
 construct complex scalar field fromr eal and imaginary parts (real-space)
 
complexScalarFieldTilde Complex (const ScalarFieldTilde &)
 convert real to complex scalar field with zero imaginary part (reciprocal-space)
 
ScalarFieldTilde O (const ScalarFieldTilde &)
 Inner product operator (diagonal in PW basis)
 
ScalarFieldTilde O (ScalarFieldTilde &&)
 Inner product operator (diagonal in PW basis)
 
complexScalarFieldTilde O (const complexScalarFieldTilde &)
 Inner product operator (diagonal in PW basis)
 
complexScalarFieldTilde O (complexScalarFieldTilde &&)
 Inner product operator (diagonal in PW basis)
 
ScalarField I (const ScalarFieldTilde &, int nThreads=0)
 Forward transform: PW basis -> real space (preserve input)
 
ScalarField I (ScalarFieldTilde &&, int nThreads=0)
 Forward transform: PW basis -> real space (destructible input)
 
complexScalarField I (const complexScalarFieldTilde &, int nThreads=0)
 Forward transform: PW basis -> real space (preserve input)
 
complexScalarField I (complexScalarFieldTilde &&, int nThreads=0)
 Forward transform: PW basis -> real space (destructible input)
 
ScalarFieldTilde J (const ScalarField &, int nThreads=0)
 Inverse transform: Real space -> PW basis.
 
complexScalarFieldTilde J (const complexScalarField &, int nThreads=0)
 Inverse transform: Real space -> PW basis (preserve input)
 
complexScalarFieldTilde J (complexScalarField &&, int nThreads=0)
 Inverse transform: Real space -> PW basis (destructible input)
 
ScalarFieldTilde Idag (const ScalarField &, int nThreads=0)
 Forward transform transpose: Real space -> PW basis.
 
complexScalarFieldTilde Idag (const complexScalarField &, int nThreads=0)
 Forward transform transpose: Real space -> PW basis (preserve input)
 
complexScalarFieldTilde Idag (complexScalarField &&, int nThreads=0)
 Forward transform transpose: Real space -> PW basis (destructible input)
 
ScalarField Jdag (const ScalarFieldTilde &, int nThreads=0)
 Inverse transform transpose: PW basis -> real space (preserve input)
 
ScalarField Jdag (ScalarFieldTilde &&, int nThreads=0)
 Inverse transform transpose: PW basis -> real space (destructible input)
 
complexScalarField Jdag (const complexScalarFieldTilde &, int nThreads=0)
 Inverse transform transpose: PW basis -> real space (preserve input)
 
complexScalarField Jdag (complexScalarFieldTilde &&, int nThreads=0)
 Inverse transform transpose: PW basis -> real space (destructible input)
 
ScalarField JdagOJ (const ScalarField &)
 Evaluate Jdag(O(J())), which avoids 2 fourier transforms in PW basis (preserve input)
 
ScalarField JdagOJ (ScalarField &&)
 Evaluate Jdag(O(J())), which avoids 2 fourier transforms in PW basis (destructible input)
 
complexScalarField JdagOJ (const complexScalarField &)
 Evaluate Jdag(O(J())), which avoids 2 fourier transforms in PW basis (preserve input)
 
complexScalarField JdagOJ (complexScalarField &&)
 Evaluate Jdag(O(J())), which avoids 2 fourier transforms in PW basis (destructible input)
 
ScalarFieldTilde L (const ScalarFieldTilde &)
 Laplacian.
 
ScalarFieldTilde L (ScalarFieldTilde &&)
 Laplacian.
 
complexScalarFieldTilde L (const complexScalarFieldTilde &)
 Laplacian.
 
complexScalarFieldTilde L (complexScalarFieldTilde &&)
 Laplacian.
 
matrix3 Lstress (const ScalarFieldTilde &X, const ScalarFieldTilde &Y)
 Lattice vector derivative of dot(X,L(Y))
 
ScalarFieldTilde Linv (const ScalarFieldTilde &)
 Inverse Laplacian.
 
ScalarFieldTilde Linv (ScalarFieldTilde &&)
 Inverse Laplacian.
 
complexScalarFieldTilde Linv (const complexScalarFieldTilde &, vector3<> *k=0)
 Inverse Laplacian.
 
complexScalarFieldTilde Linv (complexScalarFieldTilde &&, vector3<> *k=0)
 Inverse Laplacian.
 
complexScalarFieldTilde Linv (const ScalarFieldTilde &in, vector3<> *k)
 
matrix3 LinvStress (const ScalarFieldTilde &X, const ScalarFieldTilde &Y)
 Lattice vector derivative of dot(X,Linv(Y))
 
ScalarFieldTilde D (const ScalarFieldTilde &, int iDir)
 compute the gradient in the iDir'th cartesian direction
 
ScalarFieldTilde D (const ScalarFieldTilde &in, const vector3<> &dir)
 directional derivative of scalar field along (cartesian direction) dir
 
ScalarFieldTilde DD (const ScalarFieldTilde &, int iDir, int jDir)
 second derivative along iDir'th and jDir'th cartesian directions
 
void zeroNyquist (RealKernel &Gdata)
 zeros out all the nyquist components of a real G-kernel
 
void zeroNyquist (ScalarFieldTilde &Gptr)
 zeros out all the nyquist components of a G-space data array
 
void zeroNyquist (ScalarField &Rptr)
 zeros out all the nyquist components of an R-space data array
 
void removePhase (size_t N, complex *data, double &meanPhase, double &sigmaPhase, double &rmsImagErr)
 
void translate (ScalarFieldTilde &, const vector3<> &r)
 Translate by vector r in Cartesian coordinates.
 
void translate (complexScalarFieldTilde &, const vector3<> &r)
 Translate by vector r in Cartesian coordinates.
 
void multiplyBlochPhase (complexScalarField &, const vector3<> &k)
 Multiply by Block phase for wave-vector k (in reciprocal lattice coordinates)
 
ScalarField radialFunction (const GridInfo &gInfo, const RadialFunctionG &f, vector3<> r0)
 Create spherically-symmetric scalar field from radial form f, centered at lattice coordinates r0.
 
void radialFunctionG (const RadialFunctionG &f, RealKernel &Kernel)
 Create a spherically-symmetric RealKernel from its radial form f.
 
ScalarFieldTilde radialFunctionG (const GridInfo &gInfo, const RadialFunctionG &f, vector3<> r0)
 Create spherically-symmetric scalar field from radial form f, centered at lattice coordinates r0.
 
ScalarFieldTilde operator* (const RadialFunctionG &, const ScalarFieldTilde &)
 Convolve a scalar field by a radial function (preserve input)
 
ScalarFieldTilde operator* (const RadialFunctionG &, ScalarFieldTilde &&)
 Convolve a scalar field by a radial function (destructible input)
 
matrix3 convolveStress (const RadialFunctionG &w, const ScalarFieldTilde &X, const ScalarFieldTilde &Y)
 stress due to dot(X,O(w*Y))
 
ScalarField exp (const ScalarField &)
 Elementwise exponential (preserve input)
 
ScalarField exp (ScalarField &&)
 Elementwise exponential (destructible input)
 
ScalarField log (const ScalarField &)
 Elementwise logarithm (preserve input)
 
ScalarField log (ScalarField &&)
 Elementwise logarithm (destructible input)
 
ScalarField sqrt (const ScalarField &)
 Elementwise square root (preserve input)
 
ScalarField sqrt (ScalarField &&)
 Elementwise square root (destructible input)
 
ScalarField inv (const ScalarField &)
 Elementwise reciprocal (preserve input)
 
ScalarField inv (ScalarField &&)
 Elementwise reciprocal (destructible input)
 
ScalarField pow (const ScalarField &, double alpha)
 Elementwise power (preserve input)
 
ScalarField pow (ScalarField &&, double alpha)
 Elementwise power (destructible input)
 
template<class T >
Tptr clone (const Tptr &X)
 Clone (NOTE: operator= is by reference for the ScalarField classes)
 
template<class T >
Tptroperator*= (Tptr &in, double scaleFac)
 Scale.
 
template<class T >
Tptr operator* (const Tptr &in, double scaleFac)
 Scalar multiply (preserve input)
 
template<class T >
Tptr operator* (double scaleFac, const Tptr &in)
 Scalar multiply (preserve input)
 
template<class T >
Tptr operator* (Tptr &&in, double scaleFac)
 Scalar multiply (destructible input)
 
template<class T >
Tptr operator* (double scaleFac, Tptr &&in)
 Scalar multiply (destructible input)
 
template<class T >
Tptr conj (Tptr &&in)
 Generic elementwise conjugate for complex data:
 
template<class T >
Tptr conj (const Tptr &in)
 
template<class T >
Tptroperator*= (Tptr &in, const Tptr &other)
 Generic elementwise multiply for complex data:
 
ScalarFieldoperator*= (ScalarField &in, const ScalarField &other)
 Elementwise multiply for real data.
 
template<class T >
Tptr operator* (const Tptr &in1, const Tptr &in2)
 Elementwise multiply (preserve inputs)
 
template<class T >
Tptr operator* (const Tptr &in1, Tptr &&in2)
 Elementwise multiply (destructible input)
 
template<class T >
Tptr operator* (Tptr &&in1, const Tptr &in2)
 Elementwise multiply (destructible input)
 
template<class T >
Tptr operator* (Tptr &&in1, Tptr &&in2)
 Elementwise multiply (destructible inputs)
 
complexScalarFieldoperator*= (complexScalarField &, const ScalarField &)
 elementwise multiply
 
complexScalarField operator* (const complexScalarField &, const ScalarField &)
 elementwise multiply (preserve inputs)
 
complexScalarField operator* (const ScalarField &, const complexScalarField &)
 elementwise multiply (preserve inputs)
 
complexScalarField operator* (complexScalarField &&, const ScalarField &)
 elementwise multiply (destructible inputs)
 
complexScalarField operator* (const ScalarField &, complexScalarField &&)
 elementwise multiply (destructible inputs)
 
ScalarFieldTildeoperator*= (ScalarFieldTilde &, const RealKernel &)
 Elementwise multiply.
 
ScalarFieldTilde operator* (const RealKernel &, const ScalarFieldTilde &)
 Elementwise multiply (preserve inputs)
 
ScalarFieldTilde operator* (const ScalarFieldTilde &, const RealKernel &)
 Elementwise multiply (preserve inputs)
 
ScalarFieldTilde operator* (const RealKernel &, ScalarFieldTilde &&)
 Elementwise multiply (destructible input)
 
ScalarFieldTilde operator* (ScalarFieldTilde &&, const RealKernel &)
 Elementwise multiply (destructible input)
 
template<typename T >
void axpy (double alpha, const Tptr &X, Tptr &Y)
 Generic axpy for complex data types (Note: null pointers are treated as zero)
 
void axpy (double alpha, const ScalarField &X, ScalarField &Y)
 Real data Linear combine: Y += alpha * X (Note: null pointers are treated as zero)
 
template<class T >
Tptroperator+= (Tptr &in, const Tptr &other)
 Increment.
 
template<class T >
Tptroperator-= (Tptr &in, const Tptr &other)
 Decrement.
 
template<class T >
Tptr operator+ (const Tptr &in1, const Tptr &in2)
 Add (preserve inputs)
 
template<class T >
Tptr operator+ (const Tptr &in1, Tptr &&in2)
 Add (destructible input)
 
template<class T >
Tptr operator+ (Tptr &&in1, const Tptr &in2)
 Add (destructible input)
 
template<class T >
Tptr operator+ (Tptr &&in1, Tptr &&in2)
 Add (destructible inputs)
 
template<class T >
Tptr operator- (const Tptr &in1, const Tptr &in2)
 Subtract (preserve inputs)
 
template<class T >
Tptr operator- (const Tptr &in1, Tptr &&in2)
 Subtract (destructible input)
 
template<class T >
Tptr operator- (Tptr &&in1, const Tptr &in2)
 Subtract (destructible input)
 
template<class T >
Tptr operator- (Tptr &&in1, Tptr &&in2)
 Subtract (destructible inputs)
 
template<class T >
Tptr operator- (const Tptr &in)
 Negate.
 
template<class T >
Tptr operator- (Tptr &&in)
 Negate.
 
ScalarFieldoperator+= (ScalarField &, double)
 Increment by scalar.
 
ScalarField operator+ (double, const ScalarField &)
 Add scalar (preserve inputs)
 
ScalarField operator+ (const ScalarField &, double)
 Add scalar (preserve inputs)
 
ScalarField operator+ (double, ScalarField &&)
 Add scalar (destructible input)
 
ScalarField operator+ (ScalarField &&, double)
 Add scalar (destructible input)
 
ScalarFieldoperator-= (ScalarField &, double)
 Decrement by scalar.
 
ScalarField operator- (double, const ScalarField &)
 Subtract from scalar (preserve inputs)
 
ScalarField operator- (const ScalarField &, double)
 Subtract scalar (preserve inputs)
 
ScalarField operator- (double, ScalarField &&)
 Subtract from scalar (destructible input)
 
ScalarField operator- (ScalarField &&, double)
 Subtract scalar (destructible input)
 
template<typename T >
complex dot (const Tptr &X, const Tptr &Y)
 
template<typename T >
double nrm2 (const Tptr &X)
 
template<typename T >
complex sum (const Tptr &X)
 
double dot (const ScalarField &, const ScalarField &)
 Inner product.
 
double dot (const ScalarFieldTilde &, const ScalarFieldTilde &)
 Inner product.
 
double nrm2 (const ScalarField &)
 2-norm
 
double nrm2 (const ScalarFieldTilde &)
 2-norm
 
double sum (const ScalarField &)
 Sum of elements.
 
double sum (const ScalarFieldTilde &)
 Equivalent to dot() with a ScalarFieldTilde of all 1s (NOTE: sum(X) != sum(I(X)))
 
double integral (const ScalarField &)
 Integral in the unit cell (dV times sum())
 
double integral (const ScalarFieldTilde &)
 Integral in the unit cell (just fetches the G=0 component with correct prefactor)
 
complex integral (const complexScalarField &)
 Integral in the unit cell (dV times sum())
 
complex integral (const complexScalarFieldTilde &)
 Integral in the unit cell (just fetches the G=0 component with correct prefactor)
 
ScalarFieldTilde changeGrid (const ScalarFieldTilde &, const GridInfo &gInfoNew)
 
ScalarField changeGrid (const ScalarField &, const GridInfo &gInfoNew)
 
complexScalarFieldTilde changeGrid (const complexScalarFieldTilde &, const GridInfo &gInfoNew)
 
complexScalarField changeGrid (const complexScalarField &, const GridInfo &gInfoNew)
 
template<typename T >
void initZero (Tptr &X)
 
template<typename T >
void initZero (Tptr &X, const GridInfo &gInfo)
 
template<typename T >
void nullToZero (Tptr &X, const GridInfo &gInfo)
 If X is null, allocate and initialize to 0.
 
void initRandom (ScalarField &, double cap=0.0)
 initialize element-wise with a unit-normal random number (with a cap if cap>0)
 
void initRandomFlat (ScalarField &)
 initialize element-wise with a unit-flat [0:1) random number
 
void randomize (ScalarFieldTilde &)
 alternate interface needed for Mnimize
 
void initGaussianKernel (RealKernel &, double x0)
 initialize to gaussian kernel exp(-(G x0/2)^2)
 
void initTranslation (ScalarFieldTilde &, const vector3<> &r)
 initialize to translation operator exp(-i G.r)
 
ScalarFieldTilde gaussConvolve (const ScalarFieldTilde &, double sigma)
 convolve with a gaussian
 
ScalarFieldTilde gaussConvolve (ScalarFieldTilde &&, double sigma)
 convolve with a gaussian (destructible input)
 
template<typename Func , typename... Args>
void applyFuncGsq (const GridInfo &gInfo, const Func &f, Args... args)
 Evaluate a function f(i, Gsq, args...) at each point in reciprocal space indexed by i.
 
template<typename Func , typename... Args>
void applyFunc_r (const GridInfo &gInfo, const Func &f, Args... args)
 Evaluate a function f(i, r, args...) at each point in real space index by i.
 
void printStats (const ScalarField &X, const char *name, FILE *fp=stdout)
 Print mean and standard deviation of array with specified name (debug utility)
 
template<typename Callable , typename Vec >
void checkSymmetry (Callable *func, const Vec &v1, const Vec &v2, const char *funcName)
 

Detailed Description

Function Documentation

◆ checkSymmetry()

template<typename Callable , typename Vec >
void checkSymmetry ( Callable *  func,
const Vec &  v1,
const Vec &  v2,
const char *  funcName 
)

Check the symmetry of a linear operator

Template Parameters
CallableAn operator with function call signature Vec Callable(const Vec&)
VecAny operand type representing an element of a vector space

◆ dot()

template<typename T >
complex dot ( const Tptr X,
const Tptr Y 
)
Parameters
YGeneric inner product for complex types

◆ eblas_accumNorm()

void eblas_accumNorm ( int  N,
const double &  a,
const complex x,
double *  y 
)

Accumulate elementwise norm of a complex array x into y i.e. y += a x conj(x)

Parameters
NLength of array
ascale factor
xInput complex data array
yOuput real data array

◆ eblas_accumProd()

void eblas_accumProd ( int  N,
const double &  a,
const complex xU,
const complex xC,
double *  yRe,
double *  yIm 
)

Accumulate elementwise product of two complex arrays xU and xC into real and imaginary parts yRe and yIm i.e. (yRe + i yIm) += a xU conj(xC)

Parameters
NLength of array
ascale factor
xUUnconjugated input complex data array
xCConjugated input complex data array
yReOuput real-part data array
yImOuput imaginary-part data array

◆ eblas_accumProdComplex()

void eblas_accumProdComplex ( int  N,
const double &  a,
const complex xU,
const complex xC,
complex y 
)

Accumulate elementwise product of two complex arrays xU and xC into y.

Parameters
NLength of array
ascale factor
xUUnconjugated input complex data array
xCConjugated input complex data array
yOutput complex data array

◆ eblas_capMinMax()

void eblas_capMinMax ( const int  N,
double *  x,
double &  xMin,
double &  xMax,
double  capLo = -DBL_MAX,
double  capHi = +DBL_MAX 
)

Find the minimum and maximum of a data array and optionally cap it from above and/or below.

Parameters
NLength of data array
xInput data array, that might be modified if capLo or capHi are finite
xMinOn output, the minimum of the input data array
xMaxOn output, the maximum of the input data array
capLoIf finite, cap the data array on output from below at this value (no capping for the default value)
capHiIf finite, cap the data array on output from above at this value (no capping for the default value)

◆ eblas_copy()

template<typename T >
void eblas_copy ( T *  dest,
const T *  src,
int  N 
)

Copy a data array.

Template Parameters
TData type of the input and output arrays
Parameters
destDestination data pointer
srcSource data pointer
NNumber of elements of the data type T to copy

◆ eblas_div()

template<typename Ty , typename Tx >
void eblas_div ( const int  N,
const Tx *  X,
const int  incX,
Ty *  Y,
const int  incY 
)

Templated elementwise divide Y /= X for arrays X, Y.

Template Parameters
Typrimiitive data-type for array Y
Txprimiitive data-type for array X
Parameters
Nnumber of elements in X and Y
Xpointer to the first element of array X
incXstride along the X array
Ypointer to the first element of array Y
incYstride along the Y array (must be non-zero)

◆ eblas_gather_zdaxpy()

void eblas_gather_zdaxpy ( const int  Nindex,
double  a,
const int *  index,
const complex x,
complex y,
bool  conjx = false,
const complex w = 0,
bool  conjw = false 
)

Gather y += a * x(index)

Parameters
NindexLength of index array
aScale factor (real)
index0-based index array (with length at least Nindex)
xInput array that is sampled with the index array (with length at least max(index)+1)
yOutput array that is sampled consecutively (with length at least Nindex)
conjxIf true, use conj(x) instead of x
wOptional array that is sampled consecutively and multiplies x elementwise
conjwIf true, use conj(w) instead of w

◆ eblas_lincomb()

void eblas_lincomb ( const int  N,
const complex sX,
const complex X,
const int  incX,
const complex sY,
const complex Y,
const int  incY,
complex Z,
const int  incZ 
)

Elementwise linear combination Z = sX * X + sY * Y.

Parameters
NNumber of elements in X, Y and X
sXScale factor for input X
XData pointer for input X
incXPointer increment for input X
sYScale factor for input Y
YData pointer for input Y
incYPointer increment for input Y
ZData pointer for input Z
incZPointer increment for input Z

◆ eblas_mul()

template<typename Ty , typename Tx >
void eblas_mul ( const int  N,
const Tx *  X,
const int  incX,
Ty *  Y,
const int  incY 
)

Templated elementwise multiply Y *= X for arrays X, Y.

Template Parameters
Typrimiitive data-type for array Y
Txprimiitive data-type for array X
Parameters
Nnumber of elements in X and Y
Xpointer to the first element of array X
incXstride along the X array
Ypointer to the first element of array Y
incYstride along the Y array (must be non-zero)

◆ eblas_scatter_zdaxpy()

void eblas_scatter_zdaxpy ( const int  Nindex,
double  a,
const int *  index,
const complex x,
complex y,
bool  conjx = false,
const complex w = 0,
bool  conjw = false 
)

Scatter y(index) += a * x.

Parameters
NindexLength of index array
aScale factor (real)
index0-based index array (with length at least Nindex)
xInput array that is sampled consecutively (with length at least Nindex)
yOutput array that is sampled with the index array (with length at least max(index)+1)
conjxIf true, use conj(x) instead of x
wOptional array that is sampled consecutively and multiplies x elementwise
conjwIf true, use conj(w) instead of w

◆ eblas_symmetrize() [1/3]

void eblas_symmetrize ( int  N,
int  n,
const int *  symmIndex,
const int *  symmMult,
const complex phase,
complex x 
)

Symmetrize a complex array x with phase factors, using N n-fold equivalence classes in symmIndex (useful for space group symmetrization in reciprocal space)

Parameters
NLength of array x
nLength of symmetry equivalence classes
symmIndexEvery consecutive set of n indices in this array forms an equivalence class
symmMultMultiplicity per equivalence class (number of repetitions of each element in orbit)
phasePhase factors corresponding to each entry in symmIndex
xData array to be symmetrized in place

◆ eblas_symmetrize() [2/3]

void eblas_symmetrize ( int  N,
int  n,
const int *  symmIndex,
const int *  symmMult,
const complex phase,
const matrix3<> *  rotSpin,
std::vector< complex * >  x 
)

Symmetrize a quadruplet of complex arrays with phase factors, using N n-fold equivalence classes in symmIndex (useful for space group symmetrization of spin density matrices in reciprocal space)

Parameters
NLength of array x
nLength of symmetry equivalence classes
symmIndexEvery consecutive set of n indices in this array forms an equivalence class
symmMultMultiplicity per equivalence class (number of repetitions of each element in orbit)
phasePhase factors corresponding to each entry in symmIndex
rotSpinCartesian pseudo-vector (eg. spin) rotations corresponding to each symmetry operation
xData array to be symmetrized in place

◆ eblas_symmetrize() [3/3]

void eblas_symmetrize ( int  N,
int  n,
const int *  symmIndex,
double *  x 
)

Symmetrize an array x, using N n-fold equivalence classes in symmIndex.

Parameters
NLength of array x
nLength of symmetry equivalence classes
symmIndexEvery consecutive set of n indices in this array forms an equivalence class
xData array to be symmetrized in place

◆ eblas_zero()

template<typename T >
void eblas_zero ( int  N,
T *  x 
)

Zero a data array.

Template Parameters
TData type of the array
Parameters
NNumber of elements to zero
xData pointer

◆ nrm2()

template<typename T >
double nrm2 ( const Tptr X)
Parameters
XGeneric 2-norm for complex types

◆ removePhase()

void removePhase ( size_t  N,
complex data,
double &  meanPhase,
double &  sigmaPhase,
double &  rmsImagErr 
)

Convert a complex wavefunction to a real one with optimum phase choice Store the phase statistics before conversion in meanPhase and sigmaPhase and the relative rms imaginary part truncated during conversion in rmsImagErr (Useful for getting real wavefunctions in gamma point only calculations or Wannier functions)

◆ sum()

template<typename T >
complex sum ( const Tptr X)
Parameters
XGeneric sum for complex types