Solve system of equations

Syntax: | `y = GAUSSJ(m,x)`

The `GAUSSJ` function solves the system of equations `m<>x = b`, where `m` is a square matrix and `b` is a vector, and returns `x`, the vector of solutions. This function uses the Gauss-Jordan method of elimination with full pivoting. The matrix must be square. The length of the input vector must be the same as the row dimension of the matrix. The function returns a vector with the same length.

Example 1

To solve the following three equations for the three unknowns , , and :

You could use the script:

     m=[[3;10;5];[4;2;-2];[15;3;3]]
     b=[26;14;22]
     soln=gaussj(m,b)

The answer: `soln=[1.41772;-3.77215;2.4557]`

can be checked by using the outer product operator: `m<>soln` and comparing this result with the original `b` vector.

Example 2

The following script, `INVERSE.PCM`, will find the inverse of a square matrix. If you have a square matrix `M`, you could find it's inverse, `INV_M`, with the command: `@INVERSE M`

     if ( len(?1[1,*]) ~= len(?1[*,1]) ) then
       display 'input matrix must be square'
       return
     endif
     n = len(?1[1,*])
     ! identity(n) is the identity matrix of order n
     do k = [1:n]
      inv_m[1:n,k] = gaussj(?1,identity(n)[1:n,k])
     enddo