Linear Algebraic Equations

Mitra's Routines for Linear Equations

Download Linear Solver subroutines.

Sample calling program for linear solver suroutines.

Read about ALLOCATABLE statement in FORTRAN

Pseudocode for Gauss Elimination with Partial Pivoting

Gauss Elimination is for solving a system of linear algebraic equations [a]{x} = {b}. [a] is an n by n square matrix, tol is set by the user to a small number in order to detect singular matrices. If the program outputs er = -1, then the matrix is singular.

SUB Gaus ( a , b , n , x , tol , er )
er = 0
DO i = 1 , n
     s( i ) = ABS ( a ( i , 1 ) )
     DO j = 2 , n
          IF ABS ( a ( i , j ) ) > s ( i ) THEN s ( i ) = ABS ( a ( i , j ) )
     END DO
END DO
CALL Eliminate ( a , s , n , b , tol , er )
IF er .ne. -1 THEN
     CALL Substitute ( a , n , b , x )
END IF
END Gaus

SUB Eliminate ( a , s , n , b , tol , er )
DO k = 1 , n - 1
     CALL Pivot ( a , b , s , n , k )
     IF ABS ( a ( k , k ) / s ( k ) ) < tol THEN
          er = -1
          EXIT DO
     END IF
     DO i = k + 1 , n
          factor = a ( i , k ) / a ( k , k )
          DO j = k + 1 , n
               a ( i , j ) = a ( i , j ) - factor * a ( k , j )
          END DO
          b ( i ) = b ( i ) - factor * b ( k )
     END DO
END DO
IF ABS ( a ( k , k ) / s ( k ) ) < tol THEN er = -1
END Eliminate

SUB Pivot ( a , b , s , n , k )
p = k
big = ABS ( a ( k , k ) / s ( k ) )
DO ii = k + 1 , n
     dummy = ABS ( a ( k , k ) / s ( ii ) )
     IF dummy > big THEN
          big = dummy
          p = ii
     END IF
END DO
IF p .ne. k THEN
     DO jj = k , n
          dummy = a ( p , jj )
          a ( p , jj ) = a ( k , jj )
          a ( k , jj ) = dummy
     END DO
     dummy = b ( p )
     b ( p ) = b ( k )
     b ( k ) = dummy
     dummy = s ( p )
     s ( p ) = s ( k )
     s ( k ) = dummy
END IF
END Pivot

SUB Substitute ( a , b , n , x )
x( n ) = b( n ) / a ( n , n )
DO i = n - 1 , 1 , -1
     sum = 0
     DO j = i + 1 , n
          sum = sum + a ( i , j ) * x ( j )
     END DO
     x( i ) = ( b ( i ) - sum ) / a ( i , i )
END DO
END Substitute