Environment¶
This section describes the routines and data structures which help set up Elemental’s programming environment: it discusses initialization of Elemental, call stack manipulation, a custom data structure for complex data, many routines for manipulating real and complex data, a litany of custom enums, and a few useful routines for simplifying index calculations.
Build and version information¶
Every Elemental driver with proper command-line argument processing will run PrintVersion if the --version argument is used. If --build is used, then all of the below information is reported.
- void PrintVersion(std::ostream& os=std::cout )¶
Prints the Git revision, (pre-)release version, and build type. For example:
Elemental version information: Git revision: 3c6fbdaad901a554fc27a83378d63dab55af0dd3 Version: 0.81-dev Build type: PureDebug
- void PrintConfig(std::ostream& os=std::cout )¶
Prints the relevant configuration details. For example:
Elemental configuration: Math libraries: /usr/lib/liblapack.so;/usr/lib/libblas.so HAVE_F90_INTERFACE HAVE_PMRRR HAVE_MPI_REDUCE_SCATTER_BLOCK HAVE_MPI_IN_PLACE USE_BYTE_ALLGATHERS
- void PrintCCompilerInfo(std::ostream& os=std::cout )¶
Prints the relevant C compilation information. For example:
Elemental's C compiler info: CMAKE_C_COMPILER: /usr/local/bin/gcc MPI_C_COMPILER: /home/poulson/Install/bin/mpicc MPI_C_INCLUDE_PATH: /home/poulson/Install/include MPI_C_COMPILE_FLAGS: MPI_C_LINK_FLAGS: -Wl,-rpath -Wl,/home/poulson/Install/lib MPI_C_LIBRARIES: /home/poulson/Install/lib/libmpich.so;/home/poulson/Install/lib/libopa.so;/home/poulson/Install/lib/libmpl.so;/usr/lib/i386-linux-gnu/librt.so;/usr/lib/i386-linux-gnu/libpthread.so
- void PrintCxxCompilerInfo(std::ostream& os=std::cout )¶
Prints the relevant C++ compilation information. For example:
Elemental's C++ compiler info: CMAKE_CXX_COMPILER: /usr/local/bin/g++ CXX_FLAGS: -Wall MPI_CXX_COMPILER: /home/poulson/Install/bin/mpicxx MPI_CXX_INCLUDE_PATH: /home/poulson/Install/include MPI_CXX_COMPILE_FLAGS: MPI_CXX_LINK_FLAGS: -Wl,-rpath -Wl,/home/poulson/Install/lib MPI_CXX_LIBRARIES: /home/poulson/Install/lib/libmpichcxx.so;/home/poulson/Install/lib/libmpich.so;/home/poulson/Install/lib/libopa.so;/home/poulson/Install/lib/libmpl.so;/usr/lib/i386-linux-gnu/librt.so;/usr/lib/i386-linux-gnu/libpthread.so
Set up and clean up¶
- void Initialize(int& argc, char**& argv)¶
Initializes Elemental and (if necessary) MPI. The usage is very similar to MPI_Init, but the argc and argv can be directly passed in.
#include "elemental.hpp" int main( int argc, char* argv[] ) { elem::Initialize( argc, argv ); // ... elem::Finalize(); return 0; }
- void Finalize()¶
Frees all resources allocated by Elemental and (if necessary) MPI.
- bool Initialized()¶
Returns whether or not Elemental is currently initialized.
- void ReportException(std::exception& e)¶
Used for handling Elemental’s various exceptions, e.g.,
#include "elemental.hpp" int main( int argc, char* argv[] ) { elem::Initialize( argc, argv ); try { // ... } catch( std::exception& e ) { ReportException(e); } elem::Finalize(); return 0; }
Blocksize manipulation¶
- int Blocksize()¶
Return the currently chosen algorithmic blocksize. The optimal value depends on the problem size, algorithm, and architecture; the default value is 128.
- void SetBlocksize(int blocksize)¶
Change the algorithmic blocksize to the specified value.
- void PushBlocksizeStack(int blocksize)¶
It is frequently useful to temporarily change the algorithmic blocksize, so rather than having to manually store and reset the current state, one can simply push a new value onto a stack (and later pop the stack to reset the value).
- void PopBlocksizeStack()¶
Pops the stack of blocksizes. See above.
Default process grid¶
- Grid& DefaultGrid()¶
Return a process grid built over mpi::COMM_WORLD. This is typically used as a means of allowing instances of the DistMatrix<T,MC,MR> class to be constructed without having to manually specify a process grid, e.g.,
// Build a 10 x 10 distributed matrix over mpi::COMM_WORLD elem::DistMatrix<T,MC,MR> A( 10, 10 );
Call stack manipulation¶
Note
The following call stack manipulation routines are only available in non-release builds (i.e., PureDebug and HybridDebug) and are meant to allow for the call stack to be printed (via DumpCallStack()) when an exception is caught.
- void PushCallStack(std::string s)¶
Push the given routine name onto the call stack.
- void PopCallStack()¶
Remove the routine name at the top of the call stack.
- void DumpCallStack()¶
Print (and empty) the contents of the call stack.
Custom exceptions¶
- type class SingularMatrixException¶
An extension of std::runtime_error which is meant to be thrown when a singular matrix is unexpectedly encountered.
- SingularMatrixException(const char* msg="Matrix was singular")¶
Builds an instance of the exception which allows one to optionally specify the error message.
throw elem::SingularMatrixException();
- type class NonHPDMatrixException¶
An extension of std::runtime_error which is meant to be thrown when a non positive-definite Hermitian matrix is unexpectedly encountered (e.g., during Cholesky factorization).
- NonHPDMatrixException(const char* msg="Matrix was not HPD")¶
Builds an instance of the exception which allows one to optionally specify the error message.
throw elem::NonHPDMatrixException();
- type class NonHPSDMatrixException¶
An extension of std::runtime_error which is meant to be thrown when a non positive semi-definite Hermitian matrix is unexpectedly encountered (e.g., during computation of the square root of a Hermitian matrix).
- NonHPSDMatrixException(const char* msg="Matrix was not HPSD")¶
Builds an instance of the exception which allows one to optionally specify the error message.
throw elem::NonHPSDMatrixException();
Complex data¶
- type struct Complex<R>¶
- type R BaseType¶
- R real¶
The real part of the complex number
- R imag¶
The imaginary part of the complex number
- Complex()¶
This default constructor is a no-op.
- Complex(R a)¶
Construction from a real value.
- Complex(R a, R b)¶
Construction from a complex value.
- Complex(const std::complex<R>& alpha)¶
Construction from an std::complex<R> instance.
- Complex<R>& operator=(const R& alpha)¶
Assignment from a real value.
- Complex<R>& operator+=(const R& alpha)¶
Increment with a real value.
- Complex<R>& operator-=(const R& alpha)¶
Decrement with a real value.
- Complex<R>& operator*=(const R& alpha)¶
Scale with a real value.
- Complex<R>& operator/=(const R& alpha)¶
Divide with a real value.
- Complex<R>& operator=(const Complex<R>& alpha)¶
Assignment from a complex value.
- Complex<R>& operator+=(const Complex<R>& alpha)¶
Increment with a complex value.
- Complex<R>& operator-=(const Complex<R>& alpha)¶
Decrement with a complex value.
- Complex<R>& operator*=(const Complex<R>& alpha)¶
Scale with a complex value.
- Complex<R>& operator/=(const Complex<R>& alpha)¶
Divide with a complex value.
- type struct Base<F>¶
- type type¶
The underlying real datatype of the (potentially complex) datatype F. For example, typename Base<Complex<double> >::type and typename Base<double>::type are both equivalent to double. This is often extremely useful in implementing routines which are templated over real and complex datatypes but still make use of real datatypes.
- Complex<R> operator+(const Complex<R>& alpha, const Complex<R>& beta)¶
(complex,complex) addition.
- Complex<R> operator+(const Complex<R>& alpha, const R& beta)¶
(complex,real) addition.
- Complex<R> operator+(const R& alpha, const Complex<R>& beta)¶
(real,complex) addition.
- Complex<R> operator-(const Complex<R>& alpha, const Complex<R>& beta)¶
(complex,complex) subtraction.
- Complex<R> operator-(const Complex<R>& alpha, R& beta)¶
(complex,real) subtraction.
- Complex<R> operator-(const R& alpha, const Complex<R>& beta)¶
(real,complex) subtraction.
- Complex<R> operator*(const Complex<R>& alpha, const Complex<R>& beta)¶
(complex,complex) multiplication.
- Complex<R> operator*(const Complex<R>& alpha, R& beta)¶
(complex,real) multiplication.
- Complex<R> operator*(const R& alpha, const Complex<R>& beta)¶
(real,complex) multiplication.
- Complex<R> operator/(const Complex<R>& alpha, const Complex<R>& beta)¶
(complex,complex) division.
- Complex<R> operator/(const Complex<R>& alpha, const R& beta)¶
(complex,real) division.
- Complex<R> operator/(const R& alpha, const Complex<R>& beta)¶
(real,complex) division.
- Complex<R> operator+(const Complex<R>& alpha)¶
Returns alpha.
- Complex<R> operator-(const Complex<R>& alpha)¶
Returns negative alpha.
- bool operator==(const Complex<R>& alpha, const Complex<R>& beta)¶
(complex,complex) equality check.
- bool operator==(const Complex<R>& alpha, const R& beta)¶
(complex,real) equality check.
- bool operator==(const R& alpha, const Complex<R>& beta)¶
(real,complex) equality check.
- bool operator!=(const Complex<R>& alpha, const Complex<R>& beta)¶
(complex,complex) inequality check.
- bool operator!=(const Complex<R>& alpha, const R& beta)¶
(complex,real) inequality check.
- bool operator!=(const R& alpha, const Complex<R>& beta)¶
(real,complex) inequality check.
- std::ostream& operator<<(std::ostream& os, Complex<R> alpha)¶
Pretty prints alpha in the form a+bi.
- type scomplex¶
typedef Complex<float> scomplex;
- type dcomplex¶
typedef Complex<double> dcomplex;
Scalar manipulation¶
- typename Base<F>::type Abs(const F& alpha)¶
Return the absolute value of the real or complex variable \(\alpha\).
- F FastAbs(const F& alpha)¶
Return a cheaper norm of the real or complex \(\alpha\):
\[|\alpha|_{\mbox{fast}} = |\mathcal{R}(\alpha)| + |\mathcal{I}(\alpha)|\]
- F RealPart(const F& alpha)¶
Return the real part of the real or complex variable \(\alpha\).
- F ImagPart(const F& alpha)¶
Return the imaginary part of the real or complex variable \(\alpha\).
- void SetRealPart(F& alpha, typename Base<F>::type& beta)¶
Set the real part of the real or complex variable \(\alpha\) to \(\beta\).
- void SetImagPart(F& alpha, typename Base<F>::type& beta)¶
Set the imaginary part of the complex variable \(\alpha\) to \(\beta\). If \(\alpha\) has a real type, an error is thrown.
- void UpdateRealPart(F& alpha, typename Base<F>::type& beta)¶
Update the real part of the real or complex variable \(\alpha\) to \(\beta\).
- void UpdateImagPart(F& alpha, typename Base<F>::type& beta)¶
Update the imaginary part of the complex variable \(\alpha\) to \(\beta\). If \(\alpha\) has a real type, an error is thrown.
- F Conj(const F& alpha)¶
Return the complex conjugate of the real or complex variable \(\alpha\).
- F Sqrt(const F& alpha)¶
Returns the square root or the real or complex variable \(\alpha\).
- F Cos(const F& alpha)¶
Returns the cosine of the real or complex variable \(\alpha\).
- F Sin(const F& alpha)¶
Returns the sine of the real or complex variable \(\alpha\).
- F Tan(const F& alpha)¶
Returns the tangent of the real or complex variable \(\alpha\).
- F Cosh(const F& alpha)¶
Returns the hyperbolic cosine of the real or complex variable \(\alpha\).
- F Sinh(const F& alpha)¶
Returns the hyperbolic sine of the real or complex variable \(\alpha\).
- typename Base<F>::type Arg(const F& alpha)¶
Returns the argument of the real or complex variable \(\alpha\).
- Complex<R> Polar(const R& r, const R& theta=0 )¶
Returns the complex variable constructed from the polar coordinates \((r,\theta)\).
- F Exp(const F& alpha)¶
Returns the exponential of the real or complex variable \(\alpha\).
- F Pow(const F& alpha, const F& beta)¶
Returns \(\alpha^\beta\) for real or complex \(\alpha\) and \(\beta\).
- F Log(const F& alpha)¶
Returns the logarithm of the real or complex variable \(\alpha\).
Other typedefs and enums¶
- type byte¶
typedef unsigned char byte;
- type enum Conjugation¶
An enum which can be set to either CONJUGATED or UNCONJUGATED.
- type enum Distribution¶
An enum for specifying the distribution of a row or column of a distributed matrix:
- MC: Column of a standard matrix distribution
- MD: Diagonal of a standard matrix distribution
- MR: Row of a standard matrix distribution
- VC: Column-major vector distribution
- VR: Row-major vector distribution
- STAR: Redundantly stored on every process
- CIRC: Stored on a single process
- type enum ForwardOrBackward¶
An enum for specifying FORWARD or BACKWARD.
- type enum GridOrder¶
An enum for specifying either a ROW_MAJOR or COLUMN_MAJOR ordering; it is used to tune one of the algorithms in HermitianTridiag() which requires building a smaller square process grid from a rectangular process grid, as the ordering of the processes can greatly impact performance. See SetHermitianTridiagGridOrder().
- type enum LeftOrRight¶
An enum for specifying LEFT or RIGHT.
- type enum NormType¶
An enum that can be set to either
ONE_NORM:
\[\|A\|_1 = \max_{\|x\|_1=1} \|Ax\|_1 = \max_j \sum_{i=0}^{m-1} |\alpha_{i,j}|\]INFINITY_NORM:
\[\|A\|_{\infty} = \max_{\|x\|_{\infty}=1} \|Ax\|_{\infty} = \max_i \sum_{j=0}^{n-1} |\alpha_{i,j}|\]MAX_NORM:
\[\|A\|_{\mbox{max}} = \max_{i,j} |\alpha_{i,j}|\]NUCLEAR_NORM:
\[\|A\|_* = \sum_{i=0}^{\min(m,n)} \sigma_i(A)\]FROBENIUS_NORM:
\[\|A\|_F = \sqrt{\sum_{i=0}^{m-1} \sum_{j=0}^{n-1} |\alpha_{i,j}|^2} = \sum_{i=0}^{\min(m,n)} \sigma_i(A)^2\]TWO_NORM:
\[\|A\|_2 = \max_i \sigma_i(A)\]
- type enum Orientation¶
An enum for specifying whether a matrix, say \(A\), should be implicitly treated as \(A\) (NORMAL), \(A^H\) (ADJOINT), or \(A^T\) (TRANSPOSE).
- type enum UnitOrNonUnit¶
An enum for specifying either UNIT or NON_UNIT; typically used for stating whether or not a triangular matrix’s diagonal is explicitly stored (NON_UNIT) or is implicitly unit-diagonal (UNIT).
- type enum UpperOrLower¶
An enum for specifying LOWER or UPPER (triangular).
- type enum VerticalOrHorizontal¶
An enum for specifying VERTICAL or HORIZONTAL.
Indexing utilities¶
- Int Shift(Int rank, Int firstRank, Int numProcs)¶
Given a element-wise cyclic distribution over numProcs processes, where the first entry is owned by the process with rank firstRank, this routine returns the first entry owned by the process with rank rank.
- Int LocalLength(Int n, Int shift, Int numProcs)¶
Given a vector with \(n\) entries distributed over numProcs processes with shift as defined above, this routine returns the number of entries of the vector which are owned by this process.
- Int LocalLength(Int n, Int rank, Int firstRank, Int numProcs)¶
Given a vector with \(n\) entries distributed over numProcs processes, with the first entry owned by process firstRank, this routine returns the number of entries locally owned by the process with rank rank.