Direct solvers for sparse matrices, including Gaussian 
Elimination, frontal, multifrontal and supemodal, involve much 
more complicated algorithms than for dense matrices. The 
main complication is due to the need for efficiently handling 
the fill-in in the factors L and U. A typical sparse solver 
consists of four distinct steps as opposed to two in the dense 
case: 
1. An ordering step that reorders the rows and columns such 
t that the factors suffer little fill, or that the matrix has 
special structure, such as block-triangular form. 
2. An analysis step or symbolic factorization that determines 
the nonzero structures of the factors and creates suitable 
data structures for the factors. 
3. Numerical factorization that computes the L and U factors. 
4. A solve step that performs forward and back substitution 
using the factors. 
4.1 Intel PARDISO 9.0 
Intel® Math Kernel Library (Intel® MKL) offers highly 
optimized, extensively threaded math routines for scientific, 
engineering, and financial applications that require maximum 
performance. 
Multiprocessors with the PARDISO Direct Sparse Solver - an 
easy-to-use, thread-safe, high-performance, and memory- 
efficient software library licensed from the University of Basel. 
Intel MKL also includes a Conjugate Gradient iterative solver 
with a flexible reverse communication interface. 
The package PARDISO is a thread-safe, high-performance, 
robust, memory efficient and easy to use software for solving 
large sparse symmetric and unsymmetric linear systems of 
equations on shared memory multiprocessors. The solver uses a 
combination of left- and right-looking Level-3 BLAS 
supemode techniques to exploit pipelining parallelism and 
memory hierarchies. In order to improve sequential and 
parallel sparse numerical factorization performance, the 
algorithms are based on a Level-3 BLAS update and pipelining 
parallelism is exploited with a combination of left- and right 
looking supemode techniques. 
For sufficiently large problem sizes numerical experiments 
demonstrate that the scalability of the parallel algorithm is 
nearly independent of the shared-memory multiprocessing 
architecture and speedup of up to seven using eight processors 
have been observed. The approach is based on OpenMP 
directives and has been successfully tested on the following 
parallel computers: COMPAQ AlphaServer, AMD Opteron, 
SGI Origin 2000, SUN Enterprise Server, Intel Pentium 
IV/Itanium(R) with RedHat LINUX, NEC SX-5, and Cray J90 
and IBM Regatta SMPs. 
4.2 UMFPACK 
UMFPACK is a set of routines for solving unsymmetric sparse 
linear systems, Ax = b, using the Unsymmetric MultiFrontal 
method and direct sparse LU factorization. It is written in 
ANSI/ISO C, with a MATLAB interface. UMFPACK relies on 
the Level-3 Basic Linear Algebra Subprograms (dense matrix 
multiply) for its performance. This code works on Windows 
and many versions of Unix (Sun Solaris, Red Hat Linux, IBM 
AIX, SGI IRIX, and Compaq Alpha). 
UMFPACK was developed by Timothy A. Davis from the 
University of Florida, Gainesville, USA. Distributed under the 
GNU LGPL license. 
The famous math software MATLAB also use UMFPACK as 
its solver. It represent today’s highest level and many other 
solvers consider it as standard to compare. 
4.3 GRUS SPARSE SOLVER 
GSS(GRUS SPARSE SOLVER) is a software for large sparse 
matrix solving. It is written in ANSI/ISO C,can be work well 
with many sort of data and most of operating system platform. 
GSS uses new models and many new algorithms at the 
symbolic-analysis phase. GSS uses multifrontal method to 
factor A to LU, it supports complete pivoting, partial 
pivoting ,rook-pivoting and static pivoting . The tests by 
comparing the time cost of computing large mounts of matrices 
show that its average capability even better than UMFPACK in 
many field such as the max matrix order number they could 
work and their computing speed, the time GSS only cost one 
third of that UMFPACK used. 
/ >jyd xoiciou x/- an itijup: itfiflwXOK'GE3 /lda 
There is no single algorithm or software that is best for all 
types of linear systems. Some software is targeted for special 
matrices such as symmetric and positive definite, some is 
targeted for the most general cases. 
5 APPLICATION IN PHOTOGRAMMETRY 
Matrix computation are widely used in Photogrammetry such 
as Complicated coordinate translation, Adjustment 
calculation ,image enhance, image transform, image analysis 
and other image processing techniques. For a image data can 
naturally form an int matrix, image process can easily be 
translated to matrix calculation. 
Based on Meschach library the following example is core 
content of a large block indirect adjustment programme in 
photogrammetry filed to solve for the unknowns vector X of 
adjustment values in the error equations: 
V=AX-L weight P 
The coefficient matrix A is a large band sparse matrix whose 
order number may excess 10,000. L is the constant vector and 
V is error vector. A,L,P are in parameters, the others are out 
parameters. 
/Ism Li to 3V)to3vjo:;b srfr , T 2o ¿aituAoo. adi qu coisrn k v k"io 
#define MAXM 10 
sp Adj ustLS( S PM AT *A, VEC *L, SPMAT *P, VEC 
*Xout,VEC *V, double &m0) 
{ 
SPMAT *AT,*Temp,*N; 
VEC *b,*u,*ll; 
int numsteps; 
b-v_get(A->n); 
u=v_get(A->m); 
AT=sp_get( A->n, A->m,M AXM); 
for(int i=0;i<A->m;i++) 
for(int j=0;j<A->n;j-H-) 
sp set_val(AT,j ,i, sp_get_val(A,i,j));//AT=A A T 
Temp=sp_get(A->n,A->m,MAXM); 
N=sp_get( A->n, A->n,M AXM); 
sp_mltadd(AT,P,1.0,Temp); 
sp_mltadd(Temp,A,l .0,N); 
sp mv mlt(Temp,L,b); 
T emp=sp_resize(T emp, A->n, A->n); 
Temp = sp_copy(N); 
spICHfactor(T emp); 
Xout = iter_spcg(N,Temp,b,le-7,Xout,1000,&num_steps); 
//Xout out 
11 = v_get(L->dim); 
sv_mlt(-l,L,ll); 
sp_mv_mlt(A,Xout,V);//V out