# libscuff API Documentation: Utility routines

libscuff provides some simple classes
called `HMatrix`

and `HVector`

for working with numerical matrices and vectors.

This page documents a representative *subset* of the
matrix and vector functions provided; for the full API,

consult the file `include/scuff-em/libhmat.h`

that comes with the scuff-em
distribution.

# 1. Creating from scratch

Create real-valued matrices and vectors of known dimensions:

HVector *V = new HVector( 13); // create a 13-component real-valued vector HMatrix *M = new HMatrix( 5, 6); // create a 5x6 real-valued matrix

Create complex-valued matrices and vectors:

HVector *V = new HVector( 13, LHM_COMPLEX); HMatrix *M = new HMatrix( 5, 6, LHM_COMPLEX);

For symmetric or Hermitian matrices, you have the option of using
*packed storage*; this roughly halves the RAM needed to store
your matrix, but has the drawback of slowing down linear algebra
operations.

// real symmetric matrix, M_{ij} = M_{ji} HMatrix *M = new HMatrix( 5, 5, LHM_REAL, LHM_SYMMETRIC ); // complex hermitian matrix, M_{ij} = M^*_{ji} HMatrix *M = new HMatrix( 5, 5, LHM_COMPLEX, LHM_HERMITIAN ); // complex symmetric matrix, M_{ij} = M_{ji} HMatrix *M = new HMatrix( 5, 5, LHM_COMPLEX, LHM_SYMMETRIC );

# 2. Importing from text or binary files

## Text files:

Create a new vector or matrix by importing a list of numbers in an ASCII text file:

HVector *V = new HVector("MyVector.dat", LHM_TEXT ); if (V->ErrMsg) ErrExit(V->ErrMsg); HVector *M = new HMatrix("MyMatrix.dat", LHM_TEXT ); if (M->ErrMsg) ErrExit(M->ErrMsg);

Here `MyVector.dat`

should be a file with
a single number per line, while `MyMatrix.dat`

may have multiple numbers per line.
Complex numbers
are allowed. Blank lines and comments (lines beginning with a
pound sign `#`

) are skipped.
The dimension of the resulting `HVector`

,
and the number of rows of the resulting `HMatrix`

,
will be the number of non-blank non-comment lines in the file,
while the number of columns of the `HMatrix`

will
be the largest number of numbers read from any one line.

Note that, if the file import operation fails,
the constructor returns an `HVector`

or
`HMatrix`

whose `ErrMsg`

field
points to a nonempty error message. (In this case, all
other class fields in the object should be assumed to be
invalid.)

If the operation is successful, then `ErrMsg`

will be `NULL`

on return.

### HDF5 files:

Create a new vector or matrix by importing from an HDF5 binary data file:

HMatrix *M1 = new HMatrix("MyFile.hdf5", LHM_HDF5, "M1"); HMatrix *M2 = new HMatrix("MyFile.hdf5", LHM_HDF5, "M2"); HVector *V = new HVector("MyFile.hdf5", LHM_HDF5, "V"); if (M1->ErrMsg) ErrExit(M1->ErrMsg); if (M2->ErrMsg) ErrExit(M2->ErrMsg); if (V->ErrMsg) ErrExit(V->ErrMsg);

Note that the third parameter to the constructor here is the label of the dataset to within the HDF5 file.

# 3. Exporting to text or binary files

### Text files:

Write the contents of a vector or matrix to an ASCII text file:

M->ExportToText("MyMatrix.dat"); V->ExportToText("MyVector.dat");

The argument to `ExportToText`

supports `printf-`

like
semantics for inserting numbers, etc. into the file name:

int N = 3; M->ExportToText("Matrix_%i.dat", N);

### HDF5 files:

There are two calling conventions for exporting matrices and vectors to HDF5 files.

If you want each matrix and vector to be exported to a separate HDF5 file, you can say simply

M->ExportToHDF5("MyMatrix.hdf5","M"); V->ExportToHDF5("MyVector.hdf5","V");

Alternatively, you can write multiple matrices and vectors to a single HDF5 file.

void *HC = HMatrix::OpenHDF5Context("MyFile.hdf5"); M->ExportToHDF5(HC, "M"); V->ExportToHDF5(HC, "V"); HMatrix::CloseHDF5Context(HC);

# 4. Simple manipulations

** Note: All indices are zero-based.**

Get or set single entries:

HMatrix *M = new HMatrix(M, N, LHM_COMPLEX); HVector *V = new HVector(N); M->SetEntry( 3, 4, 5.6 ); M->SetEntry( 3, 4, cdouble(5.6+7.8) ); V->SetEntry( 7, 8.009 ); double D; cdouble Z; Z = M->GetEntry( 3, 4 ); D = M->GetEntryD( 3, 4 ); // discard any imaginary part Z = V->GetEntry(7); D = V->GetEntryD(8);

Get or set blocks of entries in an `HMatrix:`

HMatrix *M = new HMatrix(10, 5, LHM_COMPLEX); cdouble Col3[10]; M->GetEntriesD(":",3,Col3); cdouble M2_14[4]; M->GetEntriesD(2,"14",M2_14);

Augment individual entries:

M->AddEntry( 3, 4, 5.6 ); // the (3,4) entry gets increased by 5.6 V->AddEntry( 0, cdouble(0.0, 2.3) ); // the 0th entry gets increased by 2.3i

Replace a matrix with its conjugate or non-conjugate tranpose (no distinction for real-valued matrices):

M->Adjoint(); // conjugate transpose M->Transpose(); // non-conjugate transpose

# 5. Numerical linear algebra: LAPACK/BLAS

Multiply two matrices:

int P, Q, R; A=new HMatrix(P, Q); B=new HMatrix(Q, R); C=new HMatrix(P, R); ... A->Multiply(B, C); // set C = A*B A->Multiply(C, B, "--TransA T"); // set B = transpose(A) * C

Replace a matrix with its LU factorization:

M->LUFactorize();

Solve a single linear system using an LU-factorized matrix:

M=new HMatrix(N, N); // insert code to fill in M M->LUFactorize(); V=new HVector(N); // insert code to fill in V M->LUSolve(V); // replaces V with M^{-1} * V

Of course, we only need to call `LUFactorize`

once for
a given matrix, after which we can make any number of calls to
`LUSolve()`

with different vectors.

We can also solve multiple simultaneous systems by passing
an `HMatrix`

instead of an `HVector`

to

`LUSolve:`

M=new HMatrix(N, N); // insert code to fill in M M->LUFactorize(); R=new HMatrix(N, M); // insert code to fill in the M columns of R M->LUSolve(R); // replaces R with M^{-1} * R