Симплексный метод линейного программирования

При поиске оптимального режима симплексным методом целевая функция должна быть представлена в виде

(4.6)

а уравнения ограничений – в виде

(4.7)

x_j≥0 (j=1,2,…n). (4.8)

Принято отображать сомножители произведений, входящих в состав F_Ц, в виде векторов-строк

с=[с_1, с₂₁,…с_n] и х=[х₁, х₂,… х_n],

а коэффициенты и свободные члены в системе уравнений (4.7) – в виде векторов – столбцов

a_j

Векторы a_j называют векторами условий задачи, а вектор b называют вектором ограничений задачи. Компоненты векторов а_jв своей совокупности образуют матрицу.

A=[a₁, a₂, …a_n],

которую называют матрицей условий задачи.

В векторной форме выражение (4.6) записывается в виде скалярного произведения векторов с и х, а левая часть системы уравнений (4.7) записывается в виде произведения матрицы А на вектор х. Теперь в очень краткой форме можно записать основную задачу линейного программирования:

найти максимум функции F_Ц=сх

при условии Ах≤b, х≥0; (4.9)

или

найти минимум функции F_Ц= -сх

при условии -Ах≥-b, х≥0. (4.9')

Чтобы перейти к стандартной процедуре решения основной задачи линейного программирования симплексным методом, необходимо, чтобы среди векторов a_jбыло m единичных. Так, в системе уравнений (4.7) векторы а₁-а_m – единичные, т.е. у каждого из этих векторов имеется только одна компонента a_ij (i=1,2,…m), не равная нулю, и она равна единице. Если количество единичных векторов a_j меньше m, то в соответствующие уравнения ограничений необходимо ввести дополнительные переменные х_j, с помощью которых уравнения ограничений (4.9) будут записаны в канонической форме (4.7). В общем виде это производится методом искусственного базиса [9]. Однако в большинстве случаев недостающие переменные х_j можно просто прибавить в соответствующие уравнения ограничений, где их до того не хватало (см. пример 4.1.). Это возможно, если ограничения в условиях задачи представлены в виде неравенств.

Процедура симплексного метода основана на следующих положениях.

1. Совокупность точек х в n-мерном евклидовом пространстве Е_n, ограниченная условиями (4.9), образует область допустимых решений рассматриваемой задачи. Каждую точку х=(х₁,х₂,…х_n) из указанной области в экономических расчетах принято именовать планом (допустимым). Точка, в которой достигается либо максимум, либо минимум заданной целевой функции, называется оптимальной точкой, или оптимальным планом.

2. В пространстве Е_n область допустимых решений образуется в результате пересечения полупространств х≥0 (при n=2 – полуплоскостей) и Ах≤bили -Ах≥-b, заданных условиями (4.9), и существует в виде многогранника решений. Оптимальная точка (оптимальный план) является одной из вершин многогранника решений.

3. Любые m линейно независимых векторов a_j(j=1,2,…n, m≤n), удовлетворяющих системе уравнений (4.7), порождают точку х=(b_1,b₂,…b_m,0,…0), являющуюся вершиной многогранника решений. Такую точку принято называть угловой точкой или опорным планом. Оптимальная точка является одной из угловых точек. В этом нетрудно убедиться на примере нахождения максимального значения функции двух переменных (n=2)

(4.10)

при условиях

(4.11)

Каждое из неравенств (4.11) геометрически определяет полуплоскость в соответствии с граничными прямыми

Если система неравенств (4.11) совместна, то получим многоугольник (вырожденный многогранник, так как n=2) решений вроде приведенного на рис.4.2.

Рис. 4.2. Многоугольник решений

Область допустимых решений на рис.4.2 отмечена штриховкой ее границ. На графике рис.4.2 изображены также линии F_ц=const, соответствующие уравнению (4.10) при с₁>0 и с₂>0. Если перемещать линию F_ц в направлении вектора с=(с₁, с₂), то значение F_ц будет возрастать и достигнет максимального значения F_цm в угловой точке А. Следовательно, точка А является оптимальной, а ее координаты определяют оптимальный план решения данной задачи.

Процедура симплексного метода заключается в том, что переходят от одной угловой точки (опорного плана) к другой угловой точке (к другому опорному плану) таким образом, что значение целевой функции обязательно увеличивается (при поиске максимума целевой функции) или обязательно уменьшается (при поиске минимума). При этом указываются признаки достижения максимума или, соответственно, минимума целевой функции.

Первые m единичных векторов а_j системы уравнений (4.7)составляют базис поиска оптимального плана. Поскольку оптимальная точка является одной из угловых точек, то искать ее следует путем замены одного из векторов базиса на один из векторов а_j (j≥m+1), до того не входивших в состав базиса. После ввода в базис вектора а_j значение целевой функции можно рассчитать по формуле

(4.12)

где F_цо – предыдущее значение целевой функции;

λ – коэффициент пропорциональности, значение которого определяется в зависимости от того, какой вектор выводится из базиса.

Параметр Δ_j определяет величину и знак приращения, которое получает целевая функция после ввода в базис нового вектора а_j. Если Δ_j>0, то целевая функция убывает (F_цj <F_цо), а при Δ_j<0 целевая функция увеличивается. Знак Δ_j является критерием полезности ввода нового вектора в базис. При поиске максимума имеет смысл вводить в базис только такой вектор, которому соответствует Δ_j<0, а при поиске минимума поиск оптимального режима (оптимального плана) продолжают лишь тогда, когда находится вектор, ввод которого в базис приводит к Δ_j>0 и соответственно– к уменьшению F_ц. Значение Δ_j рассчитывается по формуле

(4.13)

где с_i-первые коэффициентов целевой функции (4.6);

а_ij – коэффициенты при х_j (j=1,2,…n) системы уравнений (4.7).

При j≤m, т.е. для m единичных векторов, входивших перед этим расчетом в состав базиса, согласно (4.13) получим Δ_j=0.

Рассмотрим более подробно поиск оптимального режима, когда последний соответствует максимуму целевой функции. Сначала перечислим признаки оптимальности режима (опорного плана).

1) Опорный план x*=(x₁^*,x₂^*,…x_m^*,0,0,…0) решения задачи (4.6)–(4.8) является оптимальным, если Δ_j ≥ 0 для любого j (j=1,2,…n).

2) Если Δ_k < 0 для некоторого j=k и среди чисел а_ik нет положительных, то целевая функция (4.6) не ограничена.

3) Если Δ_k < 0, но среди чисел a_ik есть положительные ( не все a_ik ≤ 0), то существует опорный план x΄ такой, что F_ц(x΄) > F_ц(x).

Далее составляем таблицу 4.1а. В табл. 4.1а первые m строк определяются исходными данными задачи, а показатели (m+1)-й строки вычисляют.

Таблица 4.1а

m+1 F_ц0 0 0 … 0 … 0 Δ_m+1 … Δ_k … Δ_n Исходные данные таблицы 4.1а соответствуют опорному плану х=(b₁,b₂,…b_m,0,…0).

Таблица 4.1б

m+1 F′_ц0 0 0 … Δ′_r … 0 Δ′_m+1 … 0 … Δ′_n

В (m+1)-й строке таблицы 4.1а в столбце b записывают значение F_ц0, рассчитанное по формуле (4.6), а в каждом столбце a_j (j=1, 2,…n) записывают значение Δ_j, рассчитанное по формуле (4.13). Если Δ_j<0 для некоторых индексов j и для каждого такого j по крайней мере одно из чисел a_ij положительно, то можно составить новый опорный план, при котором значение целевой функции увеличится. Чтобы перейти к улучшенному опорному плану, нужно ввести в базис вектор a_k , у которого Δ_k<0, и это Δ_k является наибольшим (по абсолютной величине) отрицательным числом среди всех Δ_j . Взамен из базиса исключается вектор a_r. Номер r соответствует элементу a_rk вектора a_k, а элемент a_rk находят по минимуму выражения b_i/a_ik для всехa_ik>0. Если же все a_ik≤0 (i=1,2,…m), то целевая функция не ограничена сверху. Число a_rk называется разрешающим элементом. Столбец a_k и строку r, на пересечении которых находится разрешающий элемент, называют направляющими (или разрешающими). Улучшенный опорный план находят методом Жордана-Гаусса. При этом положительные компоненты нового опорного плана находят по формулам:

b_i (b_rk/a_rk)∙a_ik при i ≠ r,

b′_i = (4.14)

b_rk/a_rk при i = r,

а коэффициенты разложения векторов a_j через векторы нового базиса, соответствующего новому опорному плану, - по формулам:

a_ij (a_rj/a_rk)∙a_ik при i ≠ r,

a′_ij = (4.15)

a_rj/a_rk при i = r.

После этого вместо таблицы 4.1а получится таблица 4.1б. Элементы (m+1)-й строки таблицы 4.1б находят по формулам:

F′_цо=F_цо- (b_r_k/a_rk)Δ_k , (4.16)

Δ′_j=Δ_j- (a_rj/a_rk)Δ_k , (4.17)

либо на основании их определения, см. формулы (4.6) и (4.13). Затем просматривают элементы (m+1)-й строки. Если все Δ′_j≥0, то новый опорный план оптимален. Если же имеются Δ′_j<0, то поиск оптимального плана продолжается. Заметим, что формулы (4.14)–(4.17) можно интерпретировать для визуального поиска элементов расчёта методом треугольника, когда все расчёты (кроме случая i=r) ведутся по формуле

b = a₁ – a_2∕1 a₃,

где a₁ –это либо b_i, либо a_ij ;

a_2/1 – либо b_rk/a_rk , либо a_rj/a_rk ;

a₃ – либо a_ik, либо Δ_k.

Пример 4.1. Оптимизация прокатки путем линейного программирования симплексным методом.

Определить месячную программу прокатки стали двух марок на сортовом стане, обеспечивающую максимальную прибыль.

Дано:

затраты времени по каждому сортаменту:

t₁ = 0,04 ч/т; t₂ = 0,08 ч/т;

расход энергии при прокатке:

w₁ = 8 кВт∙ч/т; w₂ = 3 кВт∙ч/т;

месячный фонд технологического времени t_м ≤ 600ч;

месячный лимит электроэнергии на текущий месяц W_м ≤ 100000 кВт∙ч.

Найти: такой вариант месячной программы прокатки, который обеспечит в текущем месяце, когда имеются ограничения по электроэнергии, максимум прибыли при условии, что прибыльность второго сортамента в 1,2 раза выше прибыльности первого сортамента стали:

F_ц = x₁ + 1,2x₂.

В выражении целевой функции F_ц прибыльность выражена в относительных единицах (о.е.). За единицу прибыли принята прибыль от 1 т проката первого сортамента.

Выпуск проката ограничивается неравенствами:

0,04x₁ + 0.08x₂ ≤ 600

8x₁ + 3x₂ ≤ 100000.

Вводим: х₃ –неиспользованное технологическое время;

х₄ –неиспользованная в пределах лимита электроэнергия;

Получим систему уравнений ограничений задачи, выраженную в канонической форме:

0,04х₁+0,08х₂+х₃ = 600;

8х₁+3х₂+х₄ = 100000.

Расчёты производим в соответствии с алгоритмом, представленным в таблицах 4.1а и 4.1б. Результаты расчётов сводим в таблицу 4.2.

Таблица 4.2

Возможный максимум прибыли составляет 13800о.е. при выпуске проката марок: х₁=12000т; х₂=1500т.