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
 

  Gamma functions
  Hermite polynomials