Systems of linear equations *Ax=b* may be divided into two classes: those with square non-degenerate *A* and those with rectangular possibly rank deficient *A*. ALGLIB package have subroutines for both types of problems.

1 Subroutines: general system

2 Subroutines: SPD/HPD system

3 Subroutines: non-square possibly singular system

4 Performance

5 Iterative refinement

6 Downloads section

The following subroutines can be used to solve real/complex systems with general form NxN non-singular matrix:

- rmatrixsolve
- rmatrixsolvem
- rmatrixlusolve
- rmatrixlusolvem
- rmatrixmixedsolve
- rmatrixmixedsolvem
- cmatrixsolve
- cmatrixsolvem
- cmatrixlusolve
- cmatrixlusolvem
- cmatrixmixedsolve
- cmatrixmixedsolvem

These subroutines can be divided into following groups:

- Those which accept real matrix (name starts with
**R**) and those which accept complex matrix (name starts with**C**). - Those which solve task with one right-hand part and those which solve task with multiple right-hand parts (name ends with
**M**). - Those which accept only
*A*, those which acceps its LU decomposition (name contains**LU**) and those which accept*both**A*and its LU decomposition (name contains**Mixed**).

The latter requires some more explanations:

- Subroutines which accept
*A*build its LU decomposition themselves in*O(N*time. Solution found may be iteratively refined (if^{ 3})**RFC**parameter is True) using extra-precision arithmetic (see below). - Subroutines which accept LU decomposition of
*A*may solve system in*O(N*time. However they can't apply iterative refinement because it requires original matrix^{ 2})*A*, not its inprecise factorization. - Subroutines which accept both
*A*and its LU decomposition combine advantages of both variants described above: high precision and high speed.

The following subroutines can be used to solve real/complex systems with symmetric/Hermitian positive definite matrix:

- spdmatrixsolve
- spdmatrixsolvem
- spdmatrixcholeskysolve
- spdmatrixcholeskysolvem
- hpdmatrixsolve
- hpdmatrixsolvem
- hpdmatrixcholeskysolve
- hpdmatrixcholeskysolvem

Naming conventions are the same as in previous section. It should be noted that subroutines from this section do not support iterative refinement.

rmatrixsolvels solves systems with rectangular possibly degenerate matrix. Is system matrix is rank deficient, solution in a least squares sense is searched for. Subroutine supports iterative refinement. Current version of ALGLIB doesn't include complex version of this subroutine or version which accept multiple right-hand parts.

Several factors influence linear solver performance:

__Performance of triangular factorization__. Highly optimized cache-oblivious algorithms are used to calculate triangular factorization. However, performance of these algorithms depends on what version of ALGLIB you are using. C++ version is the most efficient, pure NET version and versions in other languages have significantly lower performance.__Number of right-hand parts__. ALGLIB solver is optimized for systems with one right-hand part or for systems with low number of right-hand parts. If right-hand part is too large then you'll experience a performance drop becauce there is no blocking for different levels of cache for the right-hand part. In such cases you can factorize matrix yourself using cache-efficient subroutines from trfac unit, and then you can use TRSM-subroutines from ablas unit for cache-efficient solution of your system.__Iterative improvement__. Its cost is negligible in comparison with the cost of the factorization - when you solve system with one right part. However, it can require a significant part of computational resources when matrix size is very small or when you have large number of right-hand sides. In such cases you can turn it off.__Matrix copying during solution__. Algorithms which accept LU decomposition don't allocate significant amounts of additional memory. But algorithms which accept original matrix require additional memory to create its copy. However, allocation cost is negligible when matrix size is significantly less than total amount of physical RAM.

Iterative refinement is a method proposed by James H. Wilkinson to improve the accuracy of numerical solutions to systems of linear equations. After system has been solved and solution*x _{0 }* has been found, we calculate residual

**Note #1**

Iterative refinement can decrease errors which arise due to floating point round-off. It can make solution as precise as condition number of *A* allows. However, refinement result depends on accuracy with which *b-Ax _{0 }* is calculated. For iterative refinement to be worth the effort this difference must be calculated with much higher precision than that was used to find a solution. ALGLIB solver uses portable extra-precise matrix multiplier. It allows to calculate

*This article is licensed for personal use only.*

ALGLIB Project offers you two editions of ALGLIB:

**ALGLIB Free Edition**:

delivered for free

offers full set of numerical functionality

extensive algorithmic optimizations

no low level optimizations

non-commercial license

**ALGLIB Commercial Edition**:

flexible pricing

offers full set of numerical functionality

extensive algorithmic optimizations

high performance (SMP, SIMD)

commercial license with support plan

Links to download sections for Free and Commercial editions can be found below:

C++ library.

Delivered with sources.

Monolithic design.

Extreme portability.

Delivered with sources.

Monolithic design.

Extreme portability.

C# library with native kernels.

Delivered with sources.

VB.NET and IronPython wrappers.

Extreme portability.

Delivered with sources.

VB.NET and IronPython wrappers.

Extreme portability.

Delphi wrapper around C core.

Delivered as precompiled binary.

Compatible with FreePascal.

Delivered as precompiled binary.

Compatible with FreePascal.

CPython wrapper around C core.

Delivered as precompiled binary.

Delivered as precompiled binary.