I have looked long and hard for a compact Gaussian Elimination Solver that I might port to PHP and here is my favorite implementation so far:

Here is a script that illustrates usage and performs a unit test of the GaussianElimination::solve method:

<?php
/**
*  Unit test for the GaussianElimination::solve
*
* @see http://www.cs.princeton.edu/introcs/95linear/GaussianElimination.java.html
*
* Ported to PHP by Paul Meagher
* 
* @modified June 25/2007
* @modified March 25/2007
* @version 0.2
*/

// Step 1: Solve the linear equation. 

include "GaussianElimination.php";

$A = array(array(0, 1,  1),
           array(2, 4, -2),
           array(0, 3, 15));
                                 
$b = array(4, 2, 36);

$g = new GaussianElimination;

$x = $g->solve($A, $b);

// Step 2: Test whether computed and expected solution vectors are the same.

// use machine independent comparison methods
include "Utils.php";

// define expected solution vector
$expected = array(-1.0, 2.0, 2.0);

// tallys used for unit testing
$expected_success = 0;
$expected_failure = 0;

$n = count($x);

for ($i=0; $i<$n; $i++) {

  // print solution vector component
  echo $x[$i]."<br />";

  // compare computed and expected solution vectors
  if (Utils::eq($x[$i], $expected[$i]))
    $expected_success++;
  else
    $expected_failure++;  
    
}  

echo "Test of <code>GaussianElimination::solve</code> method ... "; 

// The solve method passes the test if all elements of the computed 
// solution vector match the expected solution vector

if ($expected_success == $n)
  echo "<font color='red'><b>passed</b></font>";
else 
  echo "<font color='red'><b>failed</b></font>";

?>

The output of the solve.php unit test script looks like this:

-1
2
2
Test of GaussianElimination::solve method ... passed

Adding GaussianElimination.php class to JAMA

The JAMA package does not offer a standard Gaussian Elimination solver. This is arguably because other decompositions in the package are considered superior for most of the work you might want to use a Gaussian Elimination solver on. It would nevertheless still be a good idea to have a standard Gaussian Elimination solver around for a number of reason:

• Textbooks and classrooms in linear algebra and numerical analysis still devote a lot of their curricula to discussing Gaussian Elimination so it would be good to have a standard solver to apply to text book examples and classroom projects.
• Numerical analysis code often uses a basic gaussian solver as an intermediate part of the calculation. Having a standard gaussian solver available makes porting such code easier.
• A stable guassian elimination solver might be all we need to arrive at a useful and transparent solution to a math programing problem (e.g., a spline method that uses gaussian elimination to calculate curve weights).

A natural place to put a GaussianElimination.php class is among the other JAMA decomposition methods (e.g., LUDecomposition.php, QRDecomposition, CholeskyDecompostion.php, etc...). This might be a longer term goal after the present GaussianElimination.php class has undergone further testing and development if necessary. One change that would be required would be to have the class operate on a standard JAMA matrix object ($this->A) rather than the home grown version the class currently uses. The unit testing for the class might be done a bit differently as JAMA uses vector norm tests to determine if the computed and expected solution vectors are the same.

Arbitrary precision support?

I wonder if it is necessary to add arbitrary precision support to the GaussianElimination solver? One reason for doing so is that the solver might perform better on matrices that are ill-conditioned. I would also simply be interested to see what an arbitrary precision version of GaussianElimination might look like.

Gaussian Elimination for other matrix storage formats.

The Gaussian Elimination algorithm works great with standard arrays as the storage format. If the matrix is stored in a compressed or sparse format than alternative gaussian elimination algorithms might be required that would be better adapted to the particular constraints of the matrix storage format.

Complex gaussian elimination

How would gaussian elimination on a matrix of complex numbers work? Would there be any market for gaussian elimination using matrices of complex numbers?

Visualizing and auralizing the gaussian elimination process

An interesting multimedia project might involve devising new methods to visualize and auralize the gaussian elimination process for 2D, 3D, and 4D systems of equations. One interactive mapping that might be used is the mapping between an elementrary row operation and a visual or aural quality. I have done some research and development on elementary row operations (i.e., an ElementaryRowOperations.php class) and devised a web interface to apply them to a matrix as a way to reduce the matrix to a simpler form and eventually solve it. If might interesting to explore various visualizations and auralizations of what of what each elimination step is doing to the vector space of the matrix.