Assembly and Solution - The Boundary Element Method with Programming - page 181

The procedure, which is sometimes referred to as elimination of variable n, can be
visualised as a repeated modification of the coefficients a to a * , sometimes referred to as
starring operation. We continue doing the procedure for all the equations until all the
coefficients of A below the diagonal are zero.

STEP 2: Backsubstitution
The results may now be obtained by computing the unknown from the last equation,
which involves one unknown only. The formula for computing the n-th unknown is
given by

N

1

¦

u a u F

(7.23)

n ni i n

a

nn

i n 1

The above procedure is easily converted to a subroutine. Subroutine SOLVE shown
here assumes that coefficient matrix Lhs can be stored in memory.

SUBROUTINE Solve(Lhs,Rhs,u)
!---------------------------------------------
! Solution of system of equations
! by Gauss Elimination
!---------------------------------------------
REAL(KIND=8) :: Lhs(:,:) ! Equation Left hand side
REAL(KIND=8) :: Rhs(:) ! Equation right hand side

REAL(KIND=8) :: u(:) ! Unknown

INTEGER M ! Size of system

REAL(KIND=8) :: FAC
M= UBOUND (Rhs,1)
! Reduction
Equation_n: &
DO n=1,M-1
IF (Lhs(n,n) < 1.0E-10 .AND. Lhs(n,n) > -1.0E-10) THEN
CALL Error_Message('Singular Matrix')
END IF
Equation_i: &
DO I= n + 1,M
FAC= Lhs(i,n)/Lhs(n,n)
Lhs(i, n+1 : M)= Lhs(i, n+1 : M) - Lhs(n, n+1 : M)*FAC
Rhs(i)= Rhs(i) - Rhs(n)*FAC
END DO &
Equation_i
END DO &
Equation_n
! Backsubstitution
Unknown_n: &
DO n= M,1,-1
u(n)= -1.0/Lhs(n,n)*( SUM (Lhs(n,n + 1:M)*u(n + 1:M)) - Rhs(n))
END DO Unknown_n
RETURN ; END SUBROUTINE Solve