본문 바로가기

수학적 도구/Mathematical programming

파이썬 PuLP 입문

제가 주로 하는 연구가 조합론적 최적화 쪽이다 보니 아무래도 선형 계획법(linear programming, LP)을 풀어야 할 때가 적잖이 있습니다. 그럴 때마다 예전에는 CPLEX나 Gurobi와 같은 프로그램을 C++로 불러서 해결을 했었습니다. 그렇지만 이를 사용하려면 여러 설정을 만져 줘야하는 번거로움이 있었습니다. 또, C++ 자체가 아무래도 오래된 언어다 보니 문자열 처리 등 개인적으로 불편한 부분도 많았죠.(이를 연구실 동료가 들으면 기겁하겠지만요. 하하.) 마지막으로 CPLEX나 Gurobi 모두, 비록 학술용 라이센스를 제공하기는 하지만, 상용 소프트웨어라는 점도 마음에 걸렸습니다.

 

그러던 와중에 파이썬에서 LP를 해결해 주는 오픈 소스 라이브러리가 있다는 것을 최근 알게 되었습니다. 바로 PuLP라는 라이브러리입니다.

 

Optimization with PuLP — PuLP 2.5.0 documentation

© Copyright 2009-, pulp documentation team..

coin-or.github.io

한번 써 보니 가볍고 구현하기도 쉬워서 위에서 언급한 애로사항들이 모두 해소되더군요. 물론 최적화가 필요한 실험을 해야하거나 엄청난 크기의 LP를 수행해야 한다면 아무래도 앞서 말씀드린 CPLEX나 Gurobi와 같은 상용 프로그램들이 더 좋겠지만, 가볍고 빠르게 LP를 풀어 보고 싶다면 PuLP가 제격인 것 같습니다. 여러분들에게 소개해 드릴 겸, 저도 배울 겸 해서 PuLP 입문용으로 이번 글을 작성합니다.

다운로드

거의 모든 파이썬 라이브러리가 그렇듯이 PuLP도 pip를 통해 쉽게 다운로드 받을 수 있습니다. 그냥 커맨드 창에서 아래를 입력하면 됩니다. 물론 환경 변수에 pip 경로가 지정되어 있어야 하겠죠.

 

pip install pulp

 

참고로 소스를 직접 다운 받아 설치를 할 수도 있다고 합니다. 다만 제가 해보지는 않아서 자세한 방법은 잘 모르겠네요.

기본 구조

PuLP를 사용하여 LP를 푸는 방법은 개괄적으로 아래와 같습니다.

  1. 변수(variable)를 정의한다.
  2. 문제(problem)를 정의한다.
  3. 목적 함수(objective function)를 정의한다.
  4. 제약 조건(constraints)을 정의한다.
  5. 해결한다.

예시로 다음 LP를 해결하는 PuLP 코드를 설명과 함께 단계적으로 작성해 보겠습니다.

\[ \begin{array}{lrcl} \text{minimize } & 6 x_1 + 4 x_2 + 2 x_3 & & \\ \text{subject to } & 4 x_1 + 2 x_2 + x_3 & \geq & 5, \\ & x_1 + x_2 & \geq & 3, \\ & x_2 + x_3 & \geq & 4, \\ & x_1, x_2, x_3 & \geq & 0. \end{array} \tag{1} \]

 

가장 먼저 할 일은 당연히 PuLP 라이브러리를 임포트(import)하는 일이죠. 아래 줄을 작성하시면 됩니다.

 

import pulp

 

다음은 변수를 정의해 보겠습니다. 위 LP에서 변수는 \(x_1, x_2, x_3\)가 있었습니다. 또, 모두 0보다 크거나 같다는 조건이 있었죠. 이에 해당하는 PuLP 코드는 아래와 같습니다.

 

x1 = pulp.LpVariable('x1', lowBound = 0)
x2 = pulp.LpVariable('x2', lowBound = 0)
x3 = pulp.LpVariable('x3', lowBound = 0)

 

여기서 pulp.LpVariable(…)의 첫 번째 인자(argument)는 LP를 추출할 때 나타나는 이름입니다. lowBound 인자는 해당 변수의 하한(lower bound)을 설정합니다. 반대로 상한(upper bound)을 설정하고 싶으면 upBound 인자를 정의하면 됩니다. 예를 들어 \( x_1 \)의 조건이 \( 0 \leq x_1 \leq 5 \)였다면, 아래와 같이 작성하면 됩니다.

 

x1 = pulp.LpVariable('x1', lowBound = 0, upBound = 5)

 

이제 문제를 정의하겠습니다. 여기서는 지금 풀고 싶은 문제가 최대화(maximization) 문제인지, 최소화(minimization) 문제인지를 결정합니다. 예시는 최소화 문제이므로 아래와 같이 쓰면 됩니다.

 

model = pulp.LpProblem('test_lp', pulp.LpMinimize)

 

여기서 첫 번째 인자는 똑같이 LP를 추출할 때 문제의 제목 이름에 해당합니다. 두 번째 인자는 이 문제가 최소화 문제인지 최대화 문제인지를 정해 줍니다. 만약 최대화 문제를 작성하고 싶었다면 해당 인자에 대신 pulp.LpMaximize를 넣으면 됩니다.(기본값은 최대화 문제라고 하니 정의하지 않아도 최대화 문제가 될 것입니다.)

 

이제 LP의 목적 함수를 정의하도록 하겠습니다. 예시의 목적 함수는 \(6 x_1 + 4 x_2 + 2 x_3 \)였습니다. PuLP에서 이를 추가하는 코드는 다음과 같이 매우 간단합니다.

 

model += 6 * x1 + 4 * x2 + 2 * x3

 

제약 조건도 바로 추가해 보겠습니다. 예시에서는 총 세 개의 (자명하지 않은) 제약 조건이 있었습니다. 이를 넣는 방법은 앞에서 본 목적 함수 추가 방법과 매우 유사합니다.

 

model += 4 * x1 + 2 * x2 + x3 >= 5
model += x1 + x2 >= 3
model += x2 + x3 >= 4

 

목적 함수와 제약 조건의 차이점은 부등호(혹은 등호)가 있는지 없는지로 구분됩니다. 헷갈릴 수 있으니 목적 함수를 여러 개 넣지 않도록 주의하시기 바랍니다. (목적 함수를 여러 개 넣었을 때 파레토 최적해를 구해 주는지는 실험하지 않아서 모르겠습니다.)

 

목적 함수와 제약 조건에도 이름을 붙일 수 있습니다. 바로 해당 문장 끝에 쉼표를 찍고 이름을 적으면 됩니다. 예를 들어 목적 함수의 이름을 "my_obj_func"으로 정하고 싶으면 아래와 같이 하면 됩니다.

 

model += 6 * x1 + 4 * x2 + 2 * x3, 'my_obj_func'

 

이제 문제를 풀어 보도록 하겠습니다. 이는 다음 문장을 통해 수행됩니다.

 

model.solve()

 

풀었으면 최적값(optimal value)과 최적해(optimal solution)을 얻어 내야죠. 최적값을 얻는 방법은 아래와 같습니다.

 

pulp.value(model.objective)

 

최적해, 즉 최적값을 얻기 위해 \(x_1, x_2, x_3\)이 가져야 하는 값을 얻는 방법은 다음과 같습니다.

 

x1.varValue
x2.varValue
x3.varValue

 

위 문장들은 최적값과 최적해를 출력하는 구문이 아닙니다. 해당 값을 반환하거나 담고 있을 뿐입니다. 따라서 출력을 원하시면 print(…) 함수 등을 사용하셔야 합니다.

 

이것으로 PuLP의 기본적인 구조에 대해서 알아 보았습니다. 지금까지 작성한 코드를 한 데 모으면 아래와 같아집니다.

 

import pulp

# variables
x1 = pulp.LpVariable('x1', lowBound = 0)
x2 = pulp.LpVariable('x2', lowBound = 0)
x3 = pulp.LpVariable('x3', lowBound = 0)

# model
model = pulp.LpProblem('test_lp', pulp.LpMinimize)

# objective function
model += 6 * x1 + 4 * x2 + 2 * x3

# constraints
model += 4 * x1 + 2 * x2 + x3 >= 5
model += x1 + x2 >= 3
model += x2 + x3 >= 4

# solve
model.solve()

# optimal value
print(pulp.value(model.objective))

# optimal solution
print(x1.varValue)
print(x2.varValue)
print(x3.varValue)

 

그러면 아래와 비슷한 결과를 얻을 것입니다.

 

Welcome to the CBC MILP Solver 
Version: 2.9.0
Build Date: Feb 12 2015

command line - (...)
At line 2 NAME          MODEL
At line 3 ROWS
At line 8 COLUMNS
At line 19 RHS
At line 23 BOUNDS
At line 24 ENDATA
Problem MODEL has 3 rows, 3 columns and 7 elements
Coin0008I MODEL read with 0 errors
Presolve 3 (0) rows, 3 (0) columns and 7 (0) elements
0  Obj 0 Primal inf 9.4312067 (3)
2  Obj 14
Optimal - objective value 14
Optimal objective 14 - 2 iterations time 0.002
Option for printingOptions changed from normal to all
Total time (CPU seconds):       0.00   (Wallclock seconds):       0.00
14.0
0.0
3.0
1.0

요긴한 기능들

이번 절에서는 제가 짧은 시간이나마 PuLP로 LP를 몇 개 풀어 보면서 유용하게 사용한 기능들을 소개하겠습니다.

 

LpProblem.writeLP(…)

작성한 LP를 텍스트 파일로 출력하여 줍니다. 사용하려면 LpProblem 인스턴스(instance)에서 직접 호출하면 됩니다. 즉, 위 예시에서는 아래와 같이 쓰입니다.

 

model.writeLP('out.txt')

 

이때 앞에서 변수나 제약 조건 등에 별도로 이름을 정의한 경우에는 추출된 파일에 그 이름이 그대로 적용됩니다. 다만 사용하실 때 주의할 점이 있습니다. 이름은 서로 중복되어서는 안 되며, 이름에 공백(whitespace)이 들어 가서도 안 됩니다. 공백 대신에 언더바('_')를 권장하더군요.

 

LpVariable.dicts(…)

사실 실험할 때 작성하는 LP의 변수가 한두 개인 경우는 드뭅니다. 예를 들어 어떤 그래프 문제를 푼다고 하면 정점이나 간선에 하나씩 대응하는 변수가 필요한 경우가 대부분입니다. LpVariable.dicts(…)는 이러한 변수 군(群)을 다루는 데 효과적입니다. 예를 들어, LP에 필요한 변수가 \( x_1, \cdots, x_n \)이고, 이 모든 변수들이 0보다 크거나 같아야 할 때, 이 함수를 이용해 아래와 같이 적을 수 있습니다.

 

idxlist = [i + 1 for i in range(n)]
x = pulp.LpVariable.dicts('x', idxlist, lowBound=0)

 

여기서 함수의 첫 번째 인자는 공통 이름을 나타내며, LP 추출 시 접두사(prefix)처럼 사용됩니다. 다음 인자는 인덱스 배열을 의미합니다. 함수의 이름에서 유추할 수 있듯이 파이썬의 dict와 비슷하게 숫자는 물론 다른 자료형도 키(key)가 될 수 있습니다. 이후 인자들은 LpVariable(…)과 동일합니다.

 

특정 인덱스의 변수를 호출하고 싶을 때도 파이썬의 dict와 똑같이 사용하면 됩니다. 예를 들어, \(x_1 + x_3 + x_7 \geq 6\)의 제약 조건을 만들고 싶다면 아래와 같이 쓰면 됩니다.

 

x[1] + x[3] + x[7] >= 6

 

lpSum(…)

일반적으로 LP는, 어떤 \( A \in \mathbb{R}^{m \times n} \), \( \mathbf{x}, \mathbf{c} \in \mathbb{R}^n \), \(\mathbf{b} \in \mathbb{R}^m\)에 대해 \[ \text{maximize } \mathbf{c} \cdot \mathbf{x} \text{      subject to } A \mathbf{x} \leq \mathbf{b} \tag{2} \]의 꼴로 나타내어집니다. 여기서, 목적함수는 \( \mathbf{c} \cdot \mathbf{x} = \sum_{j = 1}^n c_j x_j \)를, 제약 조건의 \(i\) 번째 행은 \( \mathbf{a}_i \cdot \mathbf{x} = \sum_{j = 1}^n a_{ij} x_j \leq b_i \)를 나타냅니다. 따라서 \( \sum_{j = 1}^n c_j x_j \)나 \( \sum_{j = 1}^n a_{ij} x_j \)를 PuLP에서 한번에 표현하는 방법이 필요합니다.

 

이를 도와주는 함수가 바로 lpSum(…)입니다. 이는 배열을 인자로 받으며 배열의 모든 원소를 덧셈으로 연결한 식을 반환합니다. 예를 들어, \( c_1, \cdots, c_n \)이 각각 \(\mathsf{c}[1], \cdots, \mathsf{c}[n]\)에 저장되어 있다고 한다면, \( \sum_{j = 1}^n c_j x_j \)는 다음과 같이 표현됩니다.

 

pulp.lpSum([c[j] * x[j] for j in idxlist])

 

같은 이치로 \( \sum_{j = 1}^n a_{ij} x_j \leq b_i \)는 다음과 같이 적을 수 있습니다.

 

pulp.lpSum([a[i][j] * x[j] for j in idxlist]) <= b[i]

 

종합하자면, 우리는 LP 2를 다음과 같이 적을 수 있습니다.

 

# variables
idxlist = [i + 1 for i in range(n)]
x = pulp.LpVariable.dicts('x', idxlist)

# model
model = pulp.LpProblem('test_lp', pulp.LpMaximize)

# objective function
model += pulp.lpSum([c[j] * x[j] for j in idxlist])

# constraints
for i in range(1, m + 1) :
    model += pulp.lpSum([a[i][j] * x[j] for j in idxlist]) <= b[i]

# solve
model. solve()

 

LpAffineExpression(…)

PuLP에서 목적 함수나 제약 조건의 각 변은 LpAffineExpression의 인스턴스입니다. 이를 활용하면 lpSum(…)처럼 한번에 만드는 것이 아니라 하나씩 단계별로 쌓아 나가면서 훨씬 섬세하게 식을 구성할 수 있습니다. 예를 들어, 아래와 같이 \[ \sum_{j = 1}^{\lfloor n / 2 \rfloor} a_{2j} x_{2j} \leq b\]

짝수에 대해서만 더한 식에 대한 제약 조건을 PuLP로 표현하면 아래와 같습니다.

 

exp = pulp.LpAffineExpression()
for j in idxlist :
    if j % 2 == 0 :
        exp += a[j] * x[j]
model += exp <= b

 

이와 동일한 작업을 앞에서 배운 lpSum(…)을 이용하여 수행할 수도 있겠지만, 개인적으로는 위 코드가 가독성도 좋고 코딩하기에도 용이하였습니다.

 

LpConstraint.pi

LP에서 쌍대성(duality)은 매우 중요한 개념으로 쌍대 문제(dual program)의 최적해는 원 문제(primal program)의 어느 제약 조건이 최적해를 얻는데 주요한 제약 조건인지를 판별하여 줍니다. 관련하여 자세한 내용은 제 이전 포스트를 참고하시기 바랍니다.

2019.07.15 - [수학적 도구/Linear programming] - [선형 계획법] 쌍대성(Duality)

 

[선형 계획법] 쌍대성(Duality)

Linear programming은 최적화 분야에서 잘 알려지고 연구되었으며 실제로도 많이 응용되는 기법입니다. 줄여서 LP라고도 하며, 우리말로는 선형계획법으로 불립니다. 이름이 뭔가 멋들어져서 괜스레

gazelle-and-cs.tistory.com

 

LP로 실험을 하다 보면, 최적해보다 오히려 중요한 제약 조건을 판단해 내는 것이 더 효과적인 통찰을 제공할 때가 있습니다. 이를 판단하기 위해서는 각 제약 조건에 대응하는 쌍대 변수(dual variable)의 값을 알아야 합니다. PuLP에서는 이를 LpConstraint.pi에 담고 있습니다. 이는 LpConstraint 클래스의 멤버 변수입니다.

 

사용 방법을 보여드리겠습니다. 기본 구조 절의 예제에서 각 제약 조건의 쌍대 변수 값을 얻는 방법은 아래와 같습니다.

 

import pulp

# variables
x1 = pulp.LpVariable('x1', lowBound = 0)
x2 = pulp.LpVariable('x2', lowBound = 0)
x3 = pulp.LpVariable('x3', lowBound = 0)

# model
model = pulp.LpProblem('test_lp', pulp.LpMinimize)

# objective function
model += 6 * x1 + 4 * x2 + 2 * x3

# constraints
c1 = 4 * x1 + 2 * x2 + x3 >= 5
c2 = x1 + x2 >= 3
c3 = x2 + x3 >= 4
model += c1
model += c2
model += c3

# solve
model.solve()

# dual values
print(c1.pi)
print(c2.pi)
print(c3.pi)

 

이를 수행해 보면 첫 번째 제약 조건은 0의 쌍대 변수 값을 갖는 반면 두 번째와 세 번째는 각각 2의 값을 갖는다는 것을 알 수 있습니다. 따라서 이 LP의 최적해를 결정하는 데 중요한 제약 조건은 두 번째와 세 번째 제약 조건이라는 사실을 유추할 수 있습니다.

마치며

이번 포스트에서는 파이썬에서 선형 계획법을 해결해 주는 오픈 소스 라이브러리인 PuLP를 사용하는 방법에 대해서 알아 보았습니다. 가볍게 LP를 구현할 수 있어서 작은 규모의 실험에서 빠르게 결과를 얻을 때 요긴하다고 생각합니다. LP를 풀어야 하는 실험이나 프로젝트를 진행하실 때 도움이 되었으면 합니다.

 

읽어 주셔서 고맙습니다. :)

'수학적 도구 > Mathematical programming' 카테고리의 다른 글

[선형 계획법] 쌍대성(Duality)  (8) 2019.07.15