BA-GMRES iterative solver with an inner iteration preconditioner
x=splspc_nrsorgmresba(A, b) x=splspc_nrsorgmresba(A, b, omg) x=splspc_nrsorgmresba(A, b, omg, maxinner) x=splspc_nrsorgmresba(A, b, omg, maxinner, maxiter) x=splspc_nrsorgmresba(A, b, omg, maxinner, maxiter, tol) [x, nbout]= splsc_nrsorgmresba(...) [x, nbout, res]= splsc_nrsorgmresba(...)
sparse matrix (size m x n) : the matrix A in the least square problem min(b-Ax)
matrix of double (size m x 1) : the vector b in the least square problem min(b-Ax)
matrix of double (size 1 x 1) : the acceleration parameter for the SOR method (default = 1.0)
matrix of integer (size 1 x 1) : the number of inner iterations (default = 1)
matrix of integer (size 1 x 1) : the maximum number of outter iterations (default = min(m,n))
matrix of double (size 1 x 1) : the threshold of tolerance for stopping the algorithm (default = 1.e-6)
matrix of double (size n x 1) : solution vector x in the least square problem min(b-Ax)
matrix of integer (size 1 x 1) : total number of outter iterations necessary to find the solution vector
matrix of double (size 1 x 1) : norm of the computed residual r=A'*b-A'*A*x
This function solves the sparse linear least squares problem min |Ax-b|, with the BA-GMRES method with NR-SOR preconditionning.
Any optional argument equal to the empty matrix [] is replaced by its default value.
If the algorithm does not converge, a warning is produced.
The GMRES method consist in solving the least square problem min(b-Ax), for A square matrix. But when A is not square, GMRES can't be directly used and need to be slightly modified.
The BA-GMRES method consist in solving the least square problem min(Bb-BAx), ie finding x such as the norm of Bb-BAx is minimal. This method is mainly useful for overdetermined problems, ie m > n for A of size m x n, because the matrix obtained by doing the product of B and A is of size n x n. However, as AB-GMRES isn't implemented for this preconditioner, underdetermined systems can be given too.
The Normal Residual Successive Over-Relaxation method (NR-SOR) is a preconditioning technique applied with an iterative method (such as GMRES or CGLS) for solving least square problems. The SOR method is an iterative technique using an acceleration parameter used to solve a linear system of equations. NR-SOR is a particular case of this method, where an updated residual is used for the computations. However, on the contrary of other preconditionning methods such as Greville or RIF, the preconditioner B isn't computed before the algorithm. Thus, in the BA-GMRES algorithm, the NR-SOR method is regulary called to roughly solve a linear system instead of using the precomputed preconditioner B. The speed and accuracy of this method strongly depend on the number of inner iterations performed at each round. The stopping criterion used is (A^T)(r_k) < tol * (A^T)(r_0).
The initial vector x0 for the iterative process has been set up to a vector of zeros.
path = fullfile(splspc_getpath(),"tests","matrices"); filename=fullfile(path,"lp_scfxm1.mtx"); A = mmread(filename); b=rand(size(A,1),1); // With default settings x=splspc_nrsorgmresba(A,b); norm(A*x-b)/norm(b) // Configure the options tol=1.e-8; maxiter=100; maxinner=4; omg=1.0; [x,nbout,res]=splspc_nrsorgmresba(A,b,omg,maxinner,maxiter,tol); // Configure only maxiter x=splspc_nrsorgmresba(A,b,[],[],maxiter); | ![]() | ![]() |
1994 - 2008 - Sparselib++ v1.7 - R. Pozo, K. Remington, and A. Lumsdaine
2011 - National Institute of Informatics - K. Hayami, K. Morikuni
2011 - DIGITEO - Michael Baudin
2011 - National Institute of Informatics - Benoit Goepfert