선형 계획 설정하기, 솔버 기반

문제를 솔버 형식으로 변환하기

이 예제에서는 솔버 기반 접근법을 사용하여 문제를 수학적 형식에서 Optimization Toolbox™ 솔버 구문으로 변환하는 방법을 보여줍니다. 이 문제는 선형 계획이지만 이 기법은 모든 솔버에 적용됩니다.

문제에 사용된 변수와 표현식은 에드거(Edgar) 및 힘멜블라우(Himmelblau)의 문헌 [1]의 예제에서 화학 공장을 가동하는 모델을 나타냅니다. 문제를 설명하는 두 개의 비디오가 있습니다.

Mathematical Modeling with Optimization, Part 1에서는 문제를 그림 형식으로 보여줍니다. 이 비디오에서는 그림에서 모델 설명에 대한 수학 표현식을 생성하는 방법을 보여줍니다.
Optimization Modeling, Part 2: Converting to Solver Form에서는 이러한 수학 표현식을 Optimization Toolbox 솔버 구문으로 변환하는 방법을 다룹니다. 이 비디오에서는 문제를 푸는 방법과 결과를 해석하는 방법을 보여줍니다.

이 예제의 나머지 부분에서는 문제를 솔버 구문으로 변환하는 방법에 대해서만 중점적으로 다룹니다. 이 예제는 비디오 Optimization Modeling, Part 2: Converting to Solver Form에서도 볼 수 있습니다. 비디오와 예제의 주요 차이점은 이 예제에서는 해시 키와 유사한 명명된 변수(또는 인덱스 변수)를 사용하는 방법을 보여준다는 점입니다. 이러한 차이점은 여러 변수를 하나의 벡터로 결합하기에 나와 있습니다.

모델 설명

비디오 Mathematical Modeling with Optimization, Part 1에서는 문제를 수학적 형식으로 변환하는 방법 중 하나는 다음이 있다고 소개합니다.

문제를 전반적으로 파악
목표 식별(최대화 또는 최소화)
(이름) 변수 식별
제약 조건 식별
제어할 수 있는 변수 확인
모든 수량을 수학 표기법으로 지정
모델이 완전하고 정확한지 확인

이 섹션에서 변수가 가지는 의미는 비디오 Mathematical Modeling with Optimization, Part 1을 참조하십시오.

최적화 문제는 목적 함수를 최소화하는 것이며, 다른 모든 표현식이 제약 조건으로 적용됩니다.

목적 함수는 다음과 같습니다.

0.002614 HPS + 0.0239 PP + 0.009825 EP.

제약 조건은 다음과 같습니다.

2500 ≤ P1 ≤ 6250
I1 ≤ 192,000
C ≤ 62,000
I1 - HE1 ≤ 132,000
I1 = LE1 + HE1 + C
1359.8 I1 = 1267.8 HE1 + 1251.4 LE1 + 192 C + 3413 P1
3000 ≤ P2 ≤ 9000
I2 ≤ 244,000
LE2 ≤ 142,000
I2 = LE2 + HE2
1359.8 I2 = 1267.8 HE2 + 1251.4 LE2 + 3413 P2
HPS = I1 + I2 + BF1
HPS = C + MPS + LPS
LPS = LE1 + LE2 + BF2
MPS = HE1 + HE2 + BF1 - BF2
P1 + P2 + PP ≥ 24,550
EP + PP ≥ 12,000
MPS ≥ 271,536
LPS ≥ 100,623
모든 변수는 양수입니다.

풀이 방법

이 최적화 문제를 풀려면 다음의 단계를 따르십시오.

솔버 선택하기
여러 변수를 하나의 벡터로 결합하기
범위 제약 조건 작성하기
선형 부등식 제약 조건 작성하기
선형 등식 제약 조건 작성하기
목적 함수 작성하기
linprog를 사용하여 문제 풀기
해 검토하기

이 단계는 비디오 Optimization Modeling, Part 2: Converting to Solver Form에서도 볼 수 있습니다.

솔버 선택하기

이 문제에 적합한 솔버를 찾으려면 최적화 의사 결정표 항목을 참조하십시오. 이 표에서는 목적 함수 유형과 제약 조건 유형을 기준으로 문제를 분류하도록 요구합니다. 이 문제에서는 목적 함수와 제약 조건이 모두 선형입니다. 의사 결정표를 따르면 linprog 솔버를 사용하면 됩니다.

Optimization Toolbox 함수가 처리하는 문제 또는 linprog 함수 도움말 페이지에서 볼 수 있듯이 linprog 솔버는 다음 형식의 문제를 풉니다.

\min_{x} f^{T} x such that {\begin{matrix} A \cdot x \leq b, \\ A e q \cdot x = b e q, \\ l b \leq x \leq u b . \end{matrix}

(1)

f^Tx는 상수 f로 구성된 행 벡터와 변수 x로 구성된 열 벡터를 곱한 값을 의미합니다. 즉, 다음과 같습니다.
f^Tx = f(1)x(1) + f(2)x(2) + ... + f(n)x(n),
여기서 n은 f의 길이입니다.
A x ≤ b는 선형 부등식을 나타냅니다. A는 k×n 행렬이고, 여기서 k는 부등식의 개수이며 n은 변수의 개수(x의 크기)입니다. b는 길이가 k인 벡터입니다. 자세한 내용은 선형 부등식 제약 조건 항목을 참조하십시오.
Aeq x = beq은 선형 등식을 나타냅니다. Aeq는 m×n 행렬이고, 여기서 m은 등식의 개수이며 n은 변수의 개수(x의 크기)입니다. beq는 길이가 m인 벡터입니다. 자세한 내용은 선형 등식 제약 조건 항목을 참조하십시오.
lb ≤ x ≤ ub는 벡터 x의 각 요소가 lb의 대응하는 요소보다 커야 하고 ub의 대응하는 요소보다 작아야 함을 의미합니다. 자세한 내용은 범위 제약 조건 항목을 참조하십시오.

함수 도움말 페이지에 표시된 것처럼 linprog 솔버의 구문은 다음과 같습니다.

[x fval] = linprog(f,A,b,Aeq,beq,lb,ub);

linprog 솔버에 대한 입력값은 수식 1에 나와 있는 행렬과 벡터입니다.

여러 변수를 하나의 벡터로 결합하기

모델 설명의 방정식에는 16개의 변수가 있습니다. 이 변수들을 하나의 벡터에 넣겠습니다. 수식 1에서 변수로 구성된 벡터의 이름은 x입니다. 순서를 결정하고 변수에서 x의 성분을 생성합니다.

다음 코드는 변수의 이름으로 구성된 셀형 배열을 사용하여 벡터를 생성합니다.

variables = {'I1','I2','HE1','HE2','LE1','LE2','C','BF1',...
    'BF2','HPS','MPS','LPS','P1','P2','PP','EP'};
N = length(variables); 
% create variables for indexing 
for v = 1:N 
   eval([variables{v},' = ', num2str(v),';']); 
end

다음 명령을 실행하면 작업 공간에 다음과 같은 명명된 변수가 생성됩니다.

Listing of variables: BF1 = 8, BF2 = 9, and so forth

이러한 명명된 변수는 x의 성분에 대한 인덱스 번호를 나타냅니다. 명명된 변수를 만들지 않아도 됩니다. 비디오 Optimization Modeling, Part 2: Converting to Solver Form에서는 x의 성분에 대한 인덱스 번호를 사용하여 간단하게 문제를 푸는 방법을 보여줍니다.

범위 제약 조건 작성하기

모델 설명의 방정식에는 하한값을 갖는 4개의 변수와 상한값을 갖는 6개의 변수가 있습니다. 하한값은 다음과 같습니다.

P1 ≥ 2500
P2 ≥ 3000
MPS ≥ 271,536
LPS ≥ 100,623.

또한, 모든 변수는 양수입니다. 이는 변수의 하한이 0임을 의미합니다.

하한 벡터 lb를 0으로 구성된 벡터로 만든 후 위의 하한 4개를 추가합니다.

lb = zeros(size(variables));
lb([P1,P2,MPS,LPS]) = ...
    [2500,3000,271536,100623];

상한이 있는 변수는 다음과 같습니다.

P1 ≤ 6250
P2 ≤ 9000
I1 ≤ 192,000
I2 ≤ 244,000
C ≤ 62,000
LE2 ≤ 142000.

상한 벡터를 Inf로 구성된 벡터로 만든 후 위의 상한 6개를 추가합니다.

ub = Inf(size(variables));
ub([P1,P2,I1,I2,C,LE2]) = ...
   [6250,9000,192000,244000,62000,142000];

선형 부등식 제약 조건 작성하기

모델 설명의 방정식에는 3개의 선형 부등식이 있습니다.

I1 - HE1 ≤ 132,000
EP + PP ≥ 12,000
P1 + P2 + PP ≥ 24,550.

A x≤b 형식의 방정식을 얻으려면 모든 변수를 부등식의 좌변에 두십시오. 위의 모든 방정식은 이미 이 형식을 가집니다. 필요한 경우, –1을 곱하여 각 부등식이 “보다 작음” 형식이 되도록 해야 합니다.

I1 - HE1 ≤ 132,000
-EP - PP ≤ -12,000
-P1 - P2 - PP ≤ -24,550.

MATLAB^® 작업 공간에서 A 행렬을 16개 변수에 대한 3개의 선형 부등식에 해당하는 3×16 영행렬로 만듭니다. 세 개의 성분을 갖는 b 벡터를 만듭니다.

A = zeros(3,16);
A(1,I1) = 1; A(1,HE1) = -1; b(1) = 132000;
A(2,EP) = -1; A(2,PP) = -1; b(2) = -12000;
A(3,[P1,P2,PP]) = [-1,-1,-1];
b(3) = -24550;

선형 등식 제약 조건 작성하기

모델 설명의 방정식에는 8개의 선형 방정식이 있습니다.

I2 = LE2 + HE2
LPS = LE1 + LE2 + BF2
HPS = I1 + I2 + BF1
HPS = C + MPS + LPS
I1 = LE1 + HE1 + C
MPS = HE1 + HE2 + BF1 - BF2
1359.8 I1 = 1267.8 HE1 + 1251.4 LE1 + 192 C + 3413 P1
1359.8 I2 = 1267.8 HE2 + 1251.4 LE2 + 3413 P2.

Aeq x=beq 형식의 방정식을 얻으려면 모든 변수를 방정식의 한쪽 변에 두십시오. 위의 방정식이 다음과 같이 됩니다.

LE2 + HE2 - I2 = 0
LE1 + LE2 + BF2 - LPS = 0
I1 + I2 + BF1 - HPS = 0
C + MPS + LPS - HPS = 0
LE1 + HE1 + C - I1 = 0
HE1 + HE2 + BF1 - BF2 - MPS = 0
1267.8 HE1 + 1251.4 LE1 + 192 C + 3413 P1 - 1359.8 I1 = 0
1267.8 HE2 + 1251.4 LE2 + 3413 P2 - 1359.8 I2 = 0.

이제 이러한 방정식에 대응되는 Aeq 행렬과 beq 벡터를 작성합니다. MATLAB 작업 공간에서 Aeq 행렬을 16개 변수에 대한 8개의 선형 방정식에 해당하는 8×16 영행렬로 만듭니다. 8개의 성분이 모두 0인 beq 벡터를 만듭니다.

Aeq = zeros(8,16); beq = zeros(8,1);
Aeq(1,[LE2,HE2,I2]) = [1,1,-1];
Aeq(2,[LE1,LE2,BF2,LPS]) = [1,1,1,-1];
Aeq(3,[I1,I2,BF1,HPS]) = [1,1,1,-1];
Aeq(4,[C,MPS,LPS,HPS]) = [1,1,1,-1];
Aeq(5,[LE1,HE1,C,I1]) = [1,1,1,-1];
Aeq(6,[HE1,HE2,BF1,BF2,MPS]) = [1,1,1,-1,-1];
Aeq(7,[HE1,LE1,C,P1,I1]) = [1267.8,1251.4,192,3413,-1359.8];
Aeq(8,[HE2,LE2,P2,I2]) = [1267.8,1251.4,3413,-1359.8];

목적 함수 작성하기

목적 함수는 다음과 같습니다.

f^Tx = 0.002614 HPS + 0.0239 PP + 0.009825 EP.

이 표현식을 x 벡터의 승수로 구성된 벡터 f로 작성합니다.

f = zeros(size(variables));
f([HPS PP EP]) = [0.002614 0.0239 0.009825];

linprog를 사용하여 문제 풀기

이제 linprog 솔버에 필요한 입력값을 얻었습니다. 다음과 같이 솔버를 호출하고 서식 지정된 형식으로 출력값을 화면에 출력하십시오.

options = optimoptions('linprog','Algorithm','dual-simplex');
[x fval] = linprog(f,A,b,Aeq,beq,lb,ub,options);
for d = 1:N
  fprintf('%12.2f \t%s\n',x(d),variables{d}) 
end
fval

다음과 같은 결과가 나타납니다.

Optimal solution found.
   136328.74 	I1
   244000.00 	I2
   128159.00 	HE1
   143377.00 	HE2
        0.00 	LE1
   100623.00 	LE2
     8169.74 	C
        0.00 	BF1
        0.00 	BF2
   380328.74 	HPS
   271536.00 	MPS
   100623.00 	LPS
     6250.00 	P1
     7060.71 	P2
    11239.29 	PP
      760.71 	EP

fval =
  1.2703e+03

해 검토하기

fval 출력값은 모든 실현가능점에서 목적 함수의 최솟값을 제공합니다.

해 벡터 x는 목적 함수가 최솟값을 가지는 점입니다. 이제 다음과 같은 동작이 일어나게 됩니다.

BF1, BF2, LE1은 0, 즉 하한입니다.
I2는 244,000, 즉 상한입니다.
f 벡터의 0이 아닌 성분은 다음과 같습니다.
- HPS — 380,328.74
- PP — 11,239.29
- EP — 760.71

비디오 Optimization Modeling, Part 2: Converting to Solver Form에서는 원래 문제 측면에서 이러한 특성에 대한 해석을 제공합니다.

참고 문헌

[1] Edgar, Thomas F., and David M. Himmelblau. Optimization of Chemical Processes. McGraw-Hill, New York, 1988.