reduced row echelon form을 구하는 Gauss-Jordan elimination 구현하기
1. reduced row echelon form
한번 쯤 다시 읽어보자
elementary row operation과 reduced row echelon form을 설명하고 있다.
elementary row operation은 1) 주어진 행렬의 i번째 행과 j번째 행을 뒤바꾸는 것
2) i번째 행에 0이 아닌 scalar를 곱하는 것
3) i번째 행에 scalar를 곱하고 j번째 행에 더해주는 것
이러한 연산을 하더라도 주어진 행렬의 rank는 변하지 않는다.
https://deepdata.tistory.com/34
reduced row echelon form은 영벡터인 행은 아래에 위치하고
각 행에서 가장 왼쪽에 있는 0이 아닌 성분은 1이며 해당 성분의 열은 1을 제외하면 모두 0이어야하고
아래에 있는 행은 위에 있는 행보다 가장 먼저 0이 아닌 성분이 오른쪽에 위치해야한다
$$\begin{pmatrix}
1 & 0 & -5\\
0 & 1 & 3\\
0 & 0 & 0\\
\end{pmatrix}$$
위 행렬이 reduced row echelon form이다.
가장 중요한 성질은 주어진 행렬의 reduced row echelon form은 유일하며, 0이 아닌 행의 수가 주어진 행렬의 rank가 된다.
2. 선형연립방정식의 해
다음과 같은 연립일차방정식을 생각해보자.
$$a_{11}x_{1} + a_{12}x_{2} + ... + a_{1m}x_{m} = b_{1}$$
$$a_{21}x_{1} + a_{22}x_{2} + ... + a_{2m}x_{m} = b_{2}$$
$$\vdots$$
$$a_{n1}x_{1} + a_{n2}x_{2} + ... + a_{nm}x_{m} = b_{n}$$
행렬로 다음과 같이 나타낸다면...
$$A = \begin{pmatrix}
a_{11} & a_{12} & ... & a_{1m} \\
a_{21} & a_{22} & ... & a_{2m} \\
\vdots & \vdots & \ddots & \vdots \\
a_{n1} & a_{n2} & ... & a_{nm} \\
\end{pmatrix}, X = \begin{pmatrix}
x_{1} \\
x_{2} \\
\vdots \\
x_{n} \end{pmatrix}, B = \begin{pmatrix}
b_{1} \\
b_{2} \\
\vdots \\
b_{n} \end{pmatrix}$$
그러면 위 연립방정식은 $$AX = B$$로 나타낼 수 있다.
방정식의 해가 유일할려면 A의 역행렬이 존재해야한다.
즉 기본적으로 A는 n*n 정방행렬이어야하며, 변수의 개수와 방정식의 개수가 일치해야한다.
3. Gauss-Jordan elimination
gaussian elimination이라고도 부른다.
연립방정식 AX = B의 해를 구하는 방법이다.
augmented matrix는 행렬 A의 맨 오른쪽 열에 벡터 B를 붙인 행렬로 다음과 같다.
이 행렬과 행렬 A의 rank가 동일하다면 주어진 선형연립방정식의 해가 존재한다는 것이 알려져있다.
반대로 rank가 서로 다르면 해가 존재하지 않는다.
또한 해가 존재하는 경우, A의 rank가 A의 열의 수 m보다 적다면 해가 무수히 많이 존재한다.
이제부터는 편의상 선형연립방정식의 해가 유일하게 존재한다는 것을 가정하자.
augmented matrix에서 elementary row operation을 적용하여 reduced row echelon form을 만들면,
마지막 열이 선형연립방정식의 유일한 해가 된다.
사실 이 과정은... 연립방정식 각각을 scalar배해서 더하고 빼는 과정으로.. 변수를 소거하는 과정과 동일하다.
2x + 3y = 5
5x + 9y = 14에서 2x+3y = 5에 3배를 하고...
6x + 9y = 15 여기서 두번째 방정식을 빼주면 x = 1이 나오고, 이를 2배 해서.. 첫번째 방정식에 빼주면
3y = 3이 되어서 y = 1이 나오듯이... 이런 과정을 수행하는거나 마찬가지다.
4. 알고리즘
구체적인 알고리즘은 다음과 같다.
1) 0행부터 n-1행까지 행렬의 행을 순회한다.
2) i번째 행의 leading element인 $a_{ii}$가 0이라면, j > i인 행에 대하여 $a_{ji}$가 0이 아닌 행과 서로 교환한다.
3) $a_{ii}$가 0이 아니라면, 각 행에 $\frac{1}{a_{ii}}$를 곱하여, leading element $a_{ii} = 1$로 만들어준다.
4) 해당 i번째 열에서 $a_{ii}$를 제외한 나머지 행 $a_{1i}, a_{2i}, ... , a_{ni}$를 모두 0으로 만들기 위하여...
i번째 행을 제외한 다른 행들을 모두 순회한다.
5) 만약 $a_{ji}$가 0이 아니라면, 이를 0으로 만들어주기 위해 i번째 행에 $a_{ji}$를 곱하고, j번째 행에 빼준다.
6) 모든 행을 순회하고 얻은 행렬의 마지막 열이 선형연립방정식의 유일한 해가 된다.
5. 예시로 이해하기
예를 들어 다음과 같은 행렬이 주어졌다.
첫번째 행부터 시작하는데 leading element가 0이므로, 다른 행과 교환하자. 두번째 행의 첫번째 원소가 2로 0이 아니므로 이와 교환한다
이제 leading element를 1로 만들기 위해 첫번째 행의 원소에 1/2를 곱해주자.
이제 1열에서 (1,1) 원소를 제외한 나머지 행의 원소들은 모두 0으로 만들어줘야한다.
세번째 행의 (3,1) 원소가 3으로 0이 아니므로, 이를 0으로 만들어줘야한다.
첫번째 행에 3배를 하고 세번째 행에 빼주면... 0으로 만들 수 있다
이제 2번째 행에 대해서도 위 과정을 반복하자.
leading element인 (2,2) 원소가 2이므로 1로 만들기 위해 행에 1/2를 곱해주고...
(1,2)와 (3,2) 원소를 모두 0으로 만들어줘야한다.
(1,2) 원소가 2이므로, 두번째 행에 2를 곱하고 첫번째 행에 빼주면 0으로 만들 수 있다.
(3,2) 원소가 -1이므로 두번째 행에 -1을 곱하고 세번째 행에 빼주면 0으로 만들 수 있다
마지막 세번째 행에서 leading element를 1로 만들어주기 위해 -2/3를 곱해주자.
역시 (1,3)과 (2,3)의 원소를 0으로 만들어줘야한다.
(1,3) 원소가 -2이므로, 3번째 행에 -2를 곱하고 첫번째 행에 빼주면.. 0으로 만들 수 있다
(2,3) 원소가 1/2이므로 3번째 행에 1/2를 곱하고 두번째 행에 빼주면 0으로 만들 수 있다
최종 결과, 마지막 열의 -16/3, -4/3, 5/3이 연립방정식의 유일한 해가 된다.
6. 연습문제
22940번: 선형 연립 방정식 (acmicpc.net)
7. 풀이
위 알고리즘을 그대로 구현하면 된다.
나눗셈을 하면서 실수연산을 필연적으로 하게 되는데,
문제 조건이 유일한 해를 가지며 정수해만을 가진다고 나와있기 때문에 마지막 결과는 정수로 출력해줘야한다.
실수연산을 하면서 오차로 인해 결과가 3.9999999999999997이나 4.000000000000000001처럼 나올 수 있는데
round()를 쓰면 가장 가까운 정수로 만들어주고 parameter를 넣지 않으면 정수형으로도 바꿔준다.
int()를 쓰면 버림을 하기 때문에 오답이 된다
또한 계수가 1이상이므로, 행을 교환할 필요는 없다.
하지만 연습을 위해 그러한 점도 구현해보자.
#선형일차연립방정식이 유일한 해를 가질때
#gaussian elimination
from sys import stdin
#i번째 행의 leading element가 0이라면, j (> i)번째 행과 교환하는 연산
def row_exchange(matrix,i):
for j in range(i+1,n):
if matrix[j][i] != 0:
matrix[i],matrix[j] = matrix[j],matrix[i]
break
return matrix
#i번째 행의 i번째 열 원소인 leading element를 1로 만들었을때,
#i번째 열의 나머지 원소를 모두 0으로 만드는 연산
def convert_to_zero(matrix,i):
for j in range(n):
if i != j:
a = matrix[j][i]
for k in range(m):
matrix[j][k] -= (matrix[i][k] * a)
return matrix
#gaussian elimination
#선형연립방정식의 유일한 정수해를 구하는 함수
def gauss_elimination(matrix):
#reduced row echelon form을 만드는 과정
for i in range(n): #모든 행을 순회하여,
if matrix[i][i] == 0: #i번째 행의 leading element가 0이라면..
matrix = row_exchange(matrix,i) #i번째 열이 0이 아닌 다른 j (> i)번째 행과 교환
#i번째 행의 leading element를 1로 만들기 위해,
#모든 0,1,..,m-1열의 원소에 1/a를 곱한다
a = matrix[i][i]
for j in range(m):
matrix[i][j] /= a
#i번째 행의 leading element가 1로 되었다면,
#i번째 열의 나머지 원소는 모두 0으로 만들어준다
matrix = convert_to_zero(matrix,i)
#reduced row echelon form을 만들었다면,
#마지막 열이 선형연립방정식의 유일한 해집합
solution = []
for i in range(n):
solution.append(round(matrix[i][-1])) #round()를 쓰면 가장 가까운 정수로 바꿔줌
return solution
n = int(stdin.readline())
matrix = []
for _ in range(n):
matrix.append(list(map(int,stdin.readline().split())))
m = len(matrix[0])
print(' '.join(map(str,gauss_elimination(matrix))))
시간복잡도는 $O(N^{3})$으로 알려져있다
Gauss-Jordan Elimination (infossm.github.io)
'알고리즘 > 선형대수학 알고리즘' 카테고리의 다른 글
행렬을 이용한 피보나치 수열 문제의 해법 (0) | 2023.03.19 |
---|---|
파이썬으로 두 행렬 곱하기 구현하는 방법 (0) | 2023.03.13 |