본 포스팅에서는 현재 코스웍으로 듣고 있는 최적화 강의 중 일부를 정리해보려고 한다. Algorithm에 대한 자세한 설명은 다루지 못할 않을 예정이지만, 전반적으로 어떻게 동작하는지 이해하는 데 도움이 되었으면 좋겠다.
Linear program (LP)를 풀기 위한 algorithm에는 크게 다음의 세 방법이 있다.
-
Kantorovich's algorithm
-
Simplex algorithm
-
Interior point method (based on strong duality)
이 중에서 simplex algorithm에 대해 알아보자.
이에 앞서, LP의 standard form을 떠올리면 다음과 같다.
min w⊺x:Ax−b≤0Cx−e=0
Simplex algorithm은 위의 식 형태와 달라 보이지만 사실은 동일한, 다음과 같은 형태를 기본으로 한다. 주목할 점은 equality constraint 대신 x의 부호에 대한 제한 조건이 추가되었다는 것이다.
max w⊺x:Ax≤bx≥0
LP의 standard form을 A, simplex algorithm은 B 로 두자. 이때 B 에서 A 로의 변환은 쉽게 받아들일 수 있다.
(B 의 두 constraint 식을 한 쪽으로 정리하면, A 의 'Ax−b≤0' inequality constraint로 한번에 정리할 수 있기 때문이다.)
그렇다면 이제 A 에서 B 로의 conversion을 알아보자. 이것이 성립하면 우리는 A 가 B 와 동치임을 알 수 있다.
min에서 max problem 으로의 변환은 단순히 부호를 flip하면 되는 문제이므로 간단하다. 여기서 핵심은 'A 에서의 equality constraint가 B 에서의 x≥0으로 어떻게 바뀌는지' 이다.
이를 간단한 예를 통해 설명해보겠다. 다음의 LP를 보자.
minx1+x2+x3:x1−x2≤0x1+2x2−3≤0x1−3x3=0
우선 두 단계를 통해 equality constraint가 없는 max problem으로 바꿔볼 것이다.
- Flipping the sign of objective function: max −(x1+x2+x3)
- Equality to inequality: '=' → '≤' + '≥'
max−x1−x2−x3:x1−x2≤0x1+2x2−3≤0x1−3x3≤03x3−x1≤0
이제 x≥0 즉, non-negativity constraint를 얻기 위한 방법을 생각해보자. 이를 위해서는 기존의 각 변수들을 새로운 변수들로 표현해주면 된다.
x1:=x4−x5(x4,x5≥0)x2:=x6−x7(x6,x7≥0)x3:=x8−x9(x8,x9≥0)
이후 x1,x2,x3를 새로운 변수 (x4∼x9)들로 표현해주면, conversion이 성립한다.
한편 simplex algorithm의 전체 과정은 복잡하기 때문에, 보다 간단한 예시로 이를 설명하려 한다.
Simplex algorithm
max w⊺x:Ax≤bx≥0
Example:
max5x1+4x2:3x1+5x2≤784x1+x2≤36x1,x2≥0
이 단계에서는 slack variable 이라는 non-negative 변수를 새로 추가하여 inequality constraint를 equality constraint로 바꿀 것이다.
즉, 위의 예시에서 slack variable (s1,s2≥0) 2개를 추가하면 문제를 다음과 같이 바꿔풀 수 있다.
maxz:z−5x1−4x2=03x1+5x2+s1=784x1+x2+s2=36x1,x2,s1,s2≥0
식을 살펴보면 우선 objective function이 z 로 대체되었고 slack variable s1,s2 가 추가되었다. 여기서 중요한 점은 78,36 자리에 위치한 숫자가 non-negative 여야 한다는 것이다. 만약 그렇지 않다면, 기존 부등식의 부호를 바꾼 다음 slack variable을 빼는 형태로 식을 구성해야 한다. (e.g. −2x1−4x2≤−34 는 2x1+4x2−s1=34 의 형태가 되도록)
Iterative algorithm (Key idea)
maxz:z−5x1−4x2=03x1+5x2+s1=784x1+x2+s2=36x1,x2,s1,s2≥0
이번에는 각 iteration 마다 다음의 단계들을 수행할 것이다.
Step 1. Initial feasible solution 을 선택
- Choose an initial feasible solution.
Step 2. Z-maximizing direction 을 따라 solution을 약간 바꿈 (perturb)
- Perturb the solution along a z-maximizing direction.
이 역시 예제를 통해 이해해보자.
Iteration #1
먼저, 다음의 두 constraint를 보자.
3x1+5x2+s1=784x1+x2+s2=36
이는 미지수 4개, 식 2개로 구성되어 있기 때문에 근의 개수가 무한히 많은 상태이다. 여기서 가장 간단한 feasible solution 을 생각해보면 x1=x2=0 일 때이고, 그때의 (x1,x2,s1,s2) = (0, 0, 78, 36) 이다. (⇒ Step 1, z=0 here)
한편 우리는 x1 이나 x2를 증가시키면 z 값이 커진다는 것을 알고 있다. 둘 (x1 or x2)을 증가시키는 방법의 수는 총 세 가지가 있는데 (하나 씩 키우는 경우 2가지 + 둘 다 키우는 경우 1가지) 이때, 각각을 키우는 2가지 case만 고려하더라도 optimal 값에 도달할 수 있다. (이는 일단 받아들이자) 이제 이 두 가지 case를 비교해보자.
- (x1,x2) = (δ,0) (δ≥0)
s1=78−3δ≥0s2=36−4δ≥0
따라서 δ≤min(26,9), δ의 최대값은 9 이므로 이때의 zmax=45 이다.
- (x1,x2) = (0,δ)
s1=78−5δ≥0s2=36−δ≥0
따라서 δ≤min(78/5,36), δ의 최대값은 78/5 이므로 이때의 zmax=62.4 이다. 즉, 이때의 zmax가 더 크기 때문에 이 perturbation 을 선택한다. (⇒ Step 2)
Iteration #2
앞선 iteration에서 얻은 두 번째 feasible solution 은 (x1,x2,s1,s2) = (0, 78/5, 0, 102/5) 이다. 이 solution을 보면 0 값을 가지고 있는 x1 이나 s1을 건들여야 할 단계임을 알 수 있는데, z 값이 x1과 x2로만 이루어져 있는 상황이라 s1이 z 값에 어떤 영향을 주는지 알기가 쉽지 않다. 따라서 이런 경우에는 Gaussian elimination을 사용한다.
+ ⇒ 5×[z−5x1−4x2=0]4×[3x1+5x2+s1=78]5z−13x1+4s1=312
이제, z-maximizing direction 을 구해보자.
maxz:5z−13x1+4s1=3123x1+5x2+s1=784x1+x2+s2=36x1,x2,s1,s2≥0
- (x1,s1) = (δ,0) (δ≥0)
5x2=78−3δ≥0x2+s2=36−4δ≥0
이때 두번째 조건은 두 개의 변수 x2,s2를 포함하고 있기 때문에 condition이 loose 할 수 있다. 따라서 현재 4x1+x2+s2의 형태로 구성된 조건을 '?⋅x1+?⋅s1+?⋅s2' 와 같이 바꿀 수 있으면 좋을 것인데, 이때도 이전과 같이 Gaussian elimination을 적용한다.
− ⇒ 5×[4x1+x2+s2=36][3x1+5x2+s1=78]17x1−s1+5s2=102
이제 아래의 식에서 Iteration #1 의 step 2 같이 perturbation 을 진행하면, (initial: (x1,x2,s1,s2) = (0, 78/5, 0, 102/5), z = 62.4)
maxz:5z−13x1+4s1=3123x1+5x2+s1=7817x1−s1+5s2=102x1,x2,s1,s2≥0
5x2=78−3δ≥05s2=102−17δ≥0
따라서 δ≤min(26,102/17), δ의 최대값은 6 이므로 이때의 zmax=78 이다. 이렇게 구한 z 값은 기존의 값 (62.4) 보다 큰 것을 확인할 수 있다.
- (x1,s1) = (0,δ)
해당 case는 s1이 증가함에 따라 z 값이 감소하기 때문에 고려하지 않는다.
Iteration #3
세 번째 feasible solution 은 (x1,x2,s1,s2) = (6, 12, 0, 0) 이다. (zmax=78)
이번에는 s1과 s2를 건들여 볼 것인데, 이전의 경우와 같이 5z−13x1+4s1 만으로는 s2가 z에 미치는 영향을 보기 어려우니 Gaussian elimination을 이용하자.
+⇒ 17×[5z−13x1+4s1=312]13×[17x1−s1+5s2=102]85z+55s1+65s2=6630
이제 기존의 문제는 아래와 같이 정리된다.
maxz:85z+55s1+65s2=66303x1+5x2+s1=7817x1−s1+5s2=102x1,x2,s1,s2≥0
그런데 이번 objective function을 보면 s1이 커지든 s2가 커지든 z는 감소한다. 따라서 iteration은 멈추고 optimal 한 (x1∗,x2∗)=(6,12), 그때의 z∗ 값은 78 이다.
Reference
카이스트 서창호 교수님 강의 - EE424: 최적화개론 lecture 6.
우연히 봤는데 좋은 글 감사합니다.