Changing the Linear Algebra Library

Kriging, which is at the root of many geostatistical algorithms requires basic linear algebra facilities, essentially matrix inversion. Most of the computing time in kriging is actually spent building and solving the kriging system, hence the importance of using an efficient linear algebra library.

The default library used by
G_{s}TL is a slightly modified version of *TNT*, the *Template Numerical Toolkit*. It is a public domain library written by Roldan Pozo, Mathematical and Computational Sciences Division, National Institute of Standards and Technology. This library was chosen because it is relatively efficient, it is easy to use and modify, and it is public domain.

However, it is easy to change the linear algebra library used by the
G_{s}TL algorithms without having to modify existing code.
The procedure is twofold. The first step is to wrap all the needed functionalities of the library into a single structure, call it `new_matrix_lib`. This structure will contain nested types (a matrix type, a vector type, ...) and static functions (for example `static void inverse(matrix& A)`), which must honor the requirements detailed in section 2.4.1.

The following is an extract from the *TNT* wrapper:

`
template<class T>
struct TNT_lib{
typedef TNT::Matrix<T> tnt_Matrix;
typedef TNT::Vector<T> tnt_Vector;
// Cholesky factorization.
static inline int cholesky(TNT::Matrix<T>& A, TNT::Matrix<T>& B){
return Cholesky_upper_factorization(A,B);
}
// LU factorization.
static inline int LU_factor(TNT::Matrix<T>& A, TNT::Vector<int>& index){
return TNT::LU_factor(A,index);
}
`
It has two nested types:

The second step is to specialize the matrix library trait class for the new library wrapper `new_matrix_lib`. The trait class is defined in <matrix_lib_traits.h>. What has been done in the case of *TNT* can serve as a model.