Next: LAPACK and BLAS. Up: Linear Systems of Algebraic Previous: Coding the LU decomposition   Contents

# Coding the forward-backward substitution

The routine accepts the LU factorization of the coeficient matrix, stored compactly in , the vector of pivots and a right hand side vector . It returns the solution of the system stored in the vector (i.e. is overwritten with the solution).

A sample substitution routine might be

```! Check if the arguments match
IF ( size(A,1).NE.size(b,1) ) THEN
PRINT*,"Error in solve: A and b must match"
RETURN
END IF
n=SIZE(A,1) ! the dimension of the matrix
!
! Permute the vector b,   b <- Pb
!
DO j=1,n-1  ! columns 1:n-1
!
IF (Pivot(j).NE.j) THEN     ! If pivot is non-diagonal
Tmp         = b(j)       ! permute elements j and pivot(j)
b(j)        = b(Pivot(j))
b(Pivot(j)) = Tmp
END IF
!
END DO

!
!  Forward Substitution,  y <- L\b
!
DO i=2,n
b(i) = b(i) - DOT_PRODUCT( A(i,1:i-1), b(1:i-1) )
END DO
!
!
!  Backward Substitution, x <- U\y
!
b(n) = b(n)/A(n,n)
DO i=n-1,1,-1                ! elements to be processed
b(i) = (b(i) - DOT_PRODUCT(A(i,i+1:n),b(i+1:n)))/A(i,i)
END DO
```