линейное программирование - это нахождение минимакса для функции f:
n
f(x₁,x₂,..,xₙ) = Σ ck * xk
k=1
в ограниченной области G пространства размерности n
ограничения:
a₁₁*x₁ + a₁₂*x₂ + .. ≤ b₁
a₂₁*x₁ + a₂₂*x₂ + .. ≤ b₂
..
набор неотрицательных переменных:
x₁ ≥ 0
x₂ ≥ 0
..
в канонической форме
c⁺ * x
A * x ≤ b
x ≥ 0
где x - вектор, который должен быть найден,
c, b - векторы заданных коэффициентов,
A - матрица заданных коэффициентов
неравенства A * x ≤ b и x ≥ 0 задают ограничения, с которыми надо оптимизировать функцию f
впервые в терминах линейного программирования в 1939г. сформулировал задачу Леонид Канторович
в 40-х годах экономист Tjalling Koopmans тоже сформулировал некоторые классические задачи экономики в терминах линейного программирования
в 1975г. Канторович и Koopmans получили премию шведского банка за работы в области экономики
в случае, если все неизвестные x должны быть целыми числами задача называется
integer programming
(IP) или integer linear programming (ILP). такие задачи считаются NP-hardесли лишь часть переменных должна быть целыми числами, то задача называется mixed integer programming (MIP). она, как правило, тоже NP-hard
особый случай - когда переменные могут принимать лишь два значения (бинарные переменные). такие переменные довольно часто используются при решении задач с прерывными условиями (discontinue conditions)
и сама функция f, и уравнения ограничений - должны быть линейными уравнениями. это означает, что не допускается перемножение переменных где-либо
правая часть уравнений (b-вектор), должен быть неизменяемым вектором:
A * x ≤ b равносильно:
A * x - b ≤ 0
если исследуемое уравнение умножить на -1, то задача максимизации станет задачей минимизации (и наоборот):
minimize 3 * x1 - x2
может также быть записана как:
maximize -3 * x1 + x2
для решения задач линейного программирования используют софт-решатель glpsol
для работы со скриптом на языке используется ключ командной строки :
glpsol --cpxlp
любая строка кода скрипта должна состять не более чем из 232 символов
начинаются с символа \ и идут до конца строки
minimize minimum min
maximize maximum max
subject to such that s.t. st. st
bounds bound
general generals gen
integer integers int
binary binaries bin
infinity inf
free
end keywords нечувствительны к регистру
могут быть :
by default all variables are non-negative, i.e. have zero lower bound and no upper bound. therefore definitions of these default bounds can be omitted in the bounds section
the integer and binary sections are intended to define some variables as integer or binary. these sections are optional and needed only in case of MIP problems. if they are specified, they must follow the bounds section or, if the latter is omitted, the constraints section
группа из ста детей, половину которых случайным образом (буквально кидая монетку для каждого), но "обязательно к исполнению", отобрали в спецпрограмму обучения, а другая половина шла в обычной программе
по факту, получили, что "показатель успеха" (мы точно знаем как определить этот фактор) для прошедших спецпрограмму 75%, а для прошедших обычную программу 45%
обозначим:
c00 лузер, без разницы какую программу он проходил
c01 добился успеха после спецпрограммы, но был бы лузером после обычной программы
c10 добился успеха после обычной программы, но был бы лузером после спецпрограммы
c11 виннер, без разницы в какой программе учился
нам необходимо определить пропорции в популяции, приводящие к максимальному успеху в целом по популяции - при заданных условиях успешности
код для решателя:
maximize success : c10 + c01 s.t. b1 : c01 + c11 = 0.75 \ это мы имеем по результатам замеров. тогда комплиментарно мы имеем также и b2 : c00 + c10 = 0.25 b3 : c10 + c11 = 0.45 \ это мы имеем по результатам замеров. тогда комплиментарно мы имеем также и b4 : c00 + c01 = 0.55 bounds 0 <= c00 <= 1 0 <= c01 <= 1 0 <= c10 <= 1 0 <= c11 <= 1 \ glpsol --cpxlp golos.lp -o output.txt
запускаемся и смотрим файл результатов:
Problem:
Rows: 4
Columns: 4
Non-zeros: 8
Status: OPTIMAL
Objective: success = 0.8 (MAXimum)
No. Row name St Activity Lower bound Upper bound Marginal
------ ------------ -- ------------- ------------- ------------- -------------
1 b1 NS 0.75 0.75 = 1
2 b2 NS 0.25 0.25 = 2
3 b3 NS 0.45 0.45 = -1
4 b4 B 0.55 0.55 =
No. Column name St Activity Lower bound Upper bound Marginal
------ ------------ -- ------------- ------------- ------------- -------------
1 c10 B 0.25 0 1
2 c01 B 0.55 0 1
3 c11 B 0.2 0 1
4 c00 NL 0 0 1 -2
Karush-Kuhn-Tucker optimality conditions:
KKT.PE: max.abs.err = 0.00e+00 on row 0
max.rel.err = 0.00e+00 on row 0
High quality
KKT.PB: max.abs.err = 0.00e+00 on row 0
max.rel.err = 0.00e+00 on row 0
High quality
KKT.DE: max.abs.err = 0.00e+00 on column 0
max.rel.err = 0.00e+00 on column 0
High quality
KKT.DB: max.abs.err = 0.00e+00 on row 0
max.rel.err = 0.00e+00 on row 0
High quality
End of output необходимо закупить два продукта в количествах, обеспечивающих минимальный уровень белков, жиров и углеводов. заданы цены продуктов. нужно найти соотношение, минимизирующее расходы :
c⁺ * x , xi > 0 , ci > 0
A * x = b , aij > 0 , bi > 0
\ жиры белки углеводы в рационе при ограничении стоимости пайка \ c цена продукта c1 = 1.1 c2 = 1.3 \ b жиры ; белки ; углеводы в сумме продуктов \ A жиры, белки, углеводы - в каждом продукте \ необходимо найти вектор x - суммарный продукт в количестве x1 x2 min f : 1.1 x1 + 1.3 x2 s.t. b1 : 1 x1 + 3 x2 > 20 b2 : 2 x1 + 2 x2 > 30 b3 : 3 x1 + 2 x2 > 40 bounds x1 > 0 x2 > 0
запускаемся и получаем в файле:
Problem:
Rows: 3
Columns: 3
Non-zeros: 9
Status: OPTIMAL
Objective: f = 22 (MINimum)
No. Row name St Activity Lower bound Upper bound Marginal
------ ------------ -- ------------- ------------- ------------- -------------
1 b1 B 40 30
2 b2 NL 40 40 0.3
3 b3 NL 50 50 0.2
No. Column name St Activity Lower bound Upper bound Marginal
------ ------------ -- ------------- ------------- ------------- -------------
1 x1 NL 0 0 0.02
2 x2 B 10 0
3 x3 B 10 0
Karush-Kuhn-Tucker optimality conditions:
KKT.PE: max.abs.err = 0.00e+00 on row 0
max.rel.err = 0.00e+00 on row 0
High quality
KKT.PB: max.abs.err = 0.00e+00 on row 0
max.rel.err = 0.00e+00 on row 0
High quality
KKT.DE: max.abs.err = 0.00e+00 on column 0
max.rel.err = 0.00e+00 on column 0
High quality
KKT.DB: max.abs.err = 0.00e+00 on row 0
max.rel.err = 0.00e+00 on row 0
High quality
End of outputкомпания по производству тары производит три вида продукции: коробки, мешки и пакеты
эта компания должна произвести за следующий календарный месяц 50 коробок и 120 мешков (есть подписанные договоры). полагаем, что пакеты находят спрос при любом произведенном количестве. по технологическим причинам компания не может произвести за месяц более 250 пакетов
у компании две фабрики: A и Б; каждая фабрика может работать не более, чем 24 дня в месяц, но и не менее чем 12 дней (КЗОТ, профсоюзы, все такое)
ежедневная производительность фабрики A составляет 2 коробки, 4 мешка и 6 пакетов
ежедневная производительность фабрики Б составляет 1 коробка, 3 мешка и 5 пакетов
полные ежедневные затраты A составляют $100
полные ежедневные затраты Б составляют $75
вы - старший экономист этой компании. сколько дней в месяц должна отработать каждая фабрика, чтобы, выполнив обязательства, минимизировать расходы?
переменные: x = дни работы фабрики A, y = дни работы фабрики Б
целевая функция: f = 100 * x + 75 * y
ограничения: обязательства и возможности
2 * x + 1 * y ≥ 50
4 * x + 3 * y ≥ 120
6 * x + 5 * y ≤ 250
12 ≤ x ≤ 24
12 ≤ y ≤ 24
кодим
\ переменные: x = дни работы фабрики A, y = дни работы фабрики Б \ целевая функция: f = 100 * x + 75 * y \ ограничения: обязательства (b1, b2) и возможности (b3) \ дни работы фабрик Minimize f : 100 x + 75 y Subject To b1 : 2 x + 1 y = 50 b2 : 4 x + 3 y = 120 b3 : 6 x + 5 y <= 250 Bounds 12 <= x <= 24 12 <= y <= 24 Integer x y
запускаемся glpsol --cpxlp factory.lp -o output.txt и смотрим файл output.txt:
Problem:
Rows: 3
Columns: 2 (2 integer, 0 binary)
Non-zeros: 6
Status: INTEGER OPTIMAL
Objective: f = 3000 (MINimum)
No. Row name Activity Lower bound Upper bound
------ ------------ ------------- ------------- -------------
1 b1 50 50 =
2 b2 120 120 =
3 b3 190 250
No. Column name Activity Lower bound Upper bound
------ ------------ ------------- ------------- -------------
1 x * 15 12 24
2 y * 20 12 24
Integer feasibility conditions:
KKT.PE: max.abs.err = 0.00e+00 on row 0
max.rel.err = 0.00e+00 on row 0
High quality
KKT.PB: max.abs.err = 0.00e+00 on row 0
max.rel.err = 0.00e+00 on row 0
High quality
End of output
\ компания должна произвести за следующий календарный месяц : 50 коробок и 120 мешков \ компания не может произвести за месяц более : 250 пакетов \ переменные: дни работы: x = фабрика A , y = фабрика Б , z = фабрика В \ целевая функция: f = 100 * x + 75 * y + 70 z \ ограничения: обязательства (b1, b2) и возможности (b3) \ ежедневная производительность фабрики A : 2 коробки, 4 мешка и 6 пакетов \ Б : 1 коробка, 3 мешка и 5 пакетов \ В : 3 коробки, 2 мешка и 4 пакета \ рабочие дни фабрик 12 - 24 Minimize f : 100 x + 75 y + 70 z Subject To b1 : 2 x + 1 y + 3 z = 50 b2 : 4 x + 3 y + 2 z = 120 b3 : 6 x + 5 y + 4 z <= 250 Bounds 12 <= x <= 24 12 <= y <= 24 12 <= z <= 24 Integer x y z
запускаемся glpsol --cpxlp factory3.lp -o output.txt и смотрим файл результатов:
Problem:
Rows: 3
Columns: 3 (3 integer, 0 binary)
Non-zeros: 9
Status: INTEGER EMPTY
Objective: f = 0 (MINimum)
No. Row name Activity Lower bound Upper bound
------ ------------ ------------- ------------- -------------
1 b1 0 50 =
2 b2 0 120 =
3 b3 0 250
No. Column name Activity Lower bound Upper bound
------ ------------ ------------- ------------- -------------
1 x * 0 12 24
2 y * 0 12 24
3 z * 0 12 24
Integer feasibility conditions:
KKT.PE: max.abs.err = 0.00e+00 on row 0
max.rel.err = 0.00e+00 on row 0
High quality
KKT.PB: max.abs.err = 1.20e+02 on row 2
max.rel.err = 9.92e-01 on row 2
SOLUTION IS INFEASIBLE
End of output вот так-то, малята... решения в целых числах - нет
попробуем решить в действительных ... уберем из файла строчки
Integer x y zкоторые в самом конце. запускаемся ...
Problem:
Rows: 3
Columns: 3
Non-zeros: 9
Status: INFEASIBLE (FINAL)
Objective: f = 1843.333333 (MINimum)
No. Row name St Activity Lower bound Upper bound Marginal
------ ------------ -- ------------- ------------- ------------- -------------
1 b1 NS 50 50 = 11.6667
2 b2 B 29.6667 120 =
3 b3 B 45.3333 250
No. Column name St Activity Lower bound Upper bound Marginal
------ ------------ -- ------------- ------------- ------------- -------------
1 x NL 12 12 24 76.6667
2 y NL 12 12 24 28.3333
3 z B -3.66667 12 24
Karush-Kuhn-Tucker optimality conditions:
KKT.PE: max.abs.err = 7.11e-15 on row 3
max.rel.err = 5.87e-17 on row 3
High quality
KKT.PB: max.abs.err = 9.03e+01 on row 2
max.rel.err = 1.21e+00 on column 3
PRIMAL SOLUTION IS INFEASIBLE
KKT.DE: max.abs.err = 3.55e-15 on column 1
max.rel.err = 7.45e-17 on column 1
High quality
KKT.DB: max.abs.err = 0.00e+00 on row 0
max.rel.err = 0.00e+00 on row 0
High quality
End of outputвона оно как! чтобы быть эффективным фабрика три должна работать минус три с половиной дня! а мужуки-то не знали!
попробуем тогда ослабить условия на равенство и "перепроизвести" коробки и мешки, в надежде найти оптимум: вернем строку с Integer, а знаки поменяем с = на =>. запускаемся...
Problem:
Rows: 3
Columns: 3 (3 integer, 0 binary)
Non-zeros: 9
Status: INTEGER OPTIMAL
Objective: f = 3150 (MINimum)
No. Row name Activity Lower bound Upper bound
------ ------------ ------------- ------------- -------------
1 b1 81 50
2 b2 123 120
3 b3 120 250
No. Column name Activity Lower bound Upper bound
------ ------------ ------------- ------------- -------------
1 x * 12 12 24
2 y * 12 12 24
3 z * 15 12 24
Integer feasibility conditions:
KKT.PE: max.abs.err = 0.00e+00 on row 0
max.rel.err = 0.00e+00 on row 0
High quality
KKT.PB: max.abs.err = 0.00e+00 on row 0
max.rel.err = 0.00e+00 on row 0
High quality
End of output
ура! а скока же коробок и мешков мы произведем в таком случае? 81 и 123. ну с мешками вроде все не так уж и страшно, а с коробками - дело швах!
ключ коммандной строки : glpsol --math
начинаются с символа # и идут до конца строки
модель задачи линейного программирования записывается на языке GMPL с использованием следующих видов объектов:
каждый объект снабжается уникальным именем, используемым для ссылок на этот объект при построении модели
объекты могут представлять собой многомерные массивы, адресуемые с помощью индексных множеств. для обращения к одному элементу такого объекта используется нотация object[index], где object - это имя объекта, а index - это перечисление значений индексных множеств. например, если demand является двумерным параметром, определенным на множестве IxJ, то конкретный элемент параметра demand будет адресоваться как demand[i,j], где i ∈ I и j ∈ J
в ряде случаев оказывается удобным описывать модель задачи и фактические данные раздельно - это позволяет использовать одну и ту же модель для решения различных экземпляров задач одного класса. с этой целью исходный файл языка GMPL состоит из двух секций: секции описания модели и секции данных. секция описания модели содержит описание всех объектов, а секция данных содержит, как правило, значения параметров для конкретного экземпляра задачи. это разделение можно усилить, поместив секции описания модели и секцию данных в различные файлы, чтобы изменение данных вообще не приводило к необходимости изменять файл, содержащий описание модели
все предложения GMPL можно разделить на две категории: декларативные и функциональные
каждое предложение должно заканчиваться точкой с запятой ;
предложение set определяет множество. если не задано пространство индексов, то set определяет простое множество, в противном случае - массив множеств по количеству элементов пространства индексов
set name alias domain, attrib, …, attrib;
name - символическое имя определяемого множества
alias - необязательная строка, задающая псевдоним множества
domain - выражение, задающее пространство индексов объекта
attrib, …, attrib - необязательные атрибуты множества:
предложение param определяет параметр задачи. если не задано пространство индексов, то param определяет один параметр (скаляр), в противном случае - массив параметров по количеству элементов пространства индексов. наиболее естественной областью применения параметров является определение числовых значений, указанных в тексте задачи (стоимости перевозок, потребности, пропускные способности) - определение параметров для этих значений позволяет записать обобщенную модель и многократно использовать ее для различных исходных данных
param name alias domain, attrib, …, attrib;
name - символическое имя определяемого параметра;
alias - необязательная строка, задающая псевдоним параметра;
domain - выражение, задающее пространство индексов параметра;
attrib, …, attrib - необязательные атрибуты параметра:
предложение var определяет переменную задачи или целый массив переменных (в случае, если задано пространство индексов). значения переменные получают в ходе решения задачи оптимизации
var name alias domain, attrib, …, attrib;
name - символическое имя определяемой переменной
alias - необязательная строка, задающая псевдоним переменной
domain –выражение, задающее пространство индексов переменной
attrib, …, attrib - необязательные атрибуты переменной:
предложение subject to определяет ограничение или целую группу ограничений (в случае, если задано пространство индексов). ключевое слово subject to может быть записано сокращенно как subj to, s.t. или даже пропущено
subject to name alias domain : expression, = expression ;
subject to name alias domain : expression, ≤ expression ;
subject to name alias domain : expression, >= expression ;
subject to name alias domain : expression, ≤ expression, ≤ expression;
subject to name alias domain : expression, >= expression, >= expression;
name - символическое имя определяемого ограничения;
alias - необязательная строка, задающая псевдоним ограничения;
domain –выражение, задающее пространство индексов ограничения;
expression - линейное выражение, используемое для вычисления компонента ограничения
предложение minimize/maximize определяет целевую функцию. при добавлении в модель нескольких целевых функций, использоваться в процессе решения задачи будет только первая из них
maximize name alias : expression ;
minimize name alias : expression ;
name - символическое имя целевой функции;
alias - необязательная строка, задающая псевдоним целевой функции;
expression - линейное выражение, используемое для вычисления значения целевой функции
для решения задачи предусмотрено необязательное предложение solve, которое может быть использовано в файле модели не более одного раза. если оно отсутствует, то предполагается, что оно находится в конце секции описания модели. в момент обработки предложения solve происходит решение задачи оптимизации, поэтому никакое из декларативных предложений не может находиться после solve. зато после обработки предложения solve переменные получают значения, соответствующие найденному решению, и к ним можно применять функциональные предложения
для вывода в GMPL предусмотрено несколько функциональных предложений. самое простое из них - display. для display может быть задано пространство индексов, а выводимые значения разделяются запятой:
display : ‘x= ‘, x, ‘y = ‘, y ;
display{o in objects} : o, take[o] ;
с помощью display можно выводить не только значения выражений, но и объекты модели (множества, ограничения, целевую функцию), однако вывод всегда производится в стандартный поток (как правило, консоль) и возможности форматирования существенно ограничены
большую гибкость предоставляет предложение printf, которое позволяет осуществлять форматированный вывод, используя синтаксис описания формата одноименной функции языка C. кроме того, printf позволяет направлять вывод в указанный файл (как с перезаписью (“>”), так и с добавлением в конец файла (“>>”)):
printf : "Hello, world!\n" ;
printf{i in 1..4} : “x[%d] = %.3f\n”, i, x[i] > “result.txt” ;
предусмотрен и высокоуровневый интерфейс для ввода-вывода табличных данных с помощью предложения table:
table name alias IN driver arg … arg: set <- [fld, …, fld] , par ~ fld, …, par ~ fld ; table name alias domain OUT driver arg … arg : expr ~ fld, … , expr ~ fld ;name –имя таблицы
описание данных расположено либо в отдельной секции файла модели, начинающейся предложением data ; и заканчивающимся предложением end ; , либо в отдельном файле (этот файл не обязательно начинать с предложения data ; )
описание данных должно задавать элементы всех множеств, определенных в секции описания модели, и значения всех указанных там параметров
задание элементов множества
set name, record, …, record ;
set name [symbol, …, symbol] , record, …, record ;
name - имя множества;
symbol, …, symbol - индексы, задающие конкретный элемент, если у множества соответствующего имени name, определено пространство индексов;
record, …, record - записи с данными:
set arcs : LED VKO KGD TJM :=
LED - + - -
VKO + - + +
KGD - - - -
TJM - - - - ; param name, record, …, record ;
param name default value, record, …, record ;
name - имя параметра;
value - необязательное значение по умолчанию для параметра;
record, …, record - записи с данными:
param trans_cost :=
[*,*,wool] :
1 2 3 :=
1 0.3 0.2 0.11
2 0.1 0.1 0.32
[*,*,steel]:
1 2 3 :=
1 0.7 0.4 0.45
2 0.2 0.2 0.64 ; пример решения первой задачи Беллмана при использовании марковских цепей. все в одном файле. если ключ -o не указан в командной строке, то файл с решением создан не будет и единственный способ увидеть решение - это использовать display, printf или table
# glpsol -m markov.mod
set States ;
set Actions ;
var value {s in States} ;
param transition_probability {s in States , a in Actions , s2 in States} ;
param reward {s in States , a in Actions} ;
param discount_factor ;
minimize
total: sum {s in States} value[s] ;
s.t.
bellman {s in States , a in Actions}:
value[s] >=
reward[s,a] +
sum {s2 in States} discount_factor *
transition_probability[s , a , s2] *
value[s2] ;
solve ;
display : total ;
end ;
data ;
set States := good deteriorating broken ;
set Actions := maintain ignore ;
param transition_probability:=
[* , * , good]:
maintain ignore :=
good 1 .5
deteriorating .9 0
broken .2 0
[* , * , deteriorating]:
maintain ignore :=
good 0 .5
deteriorating .1 .5
broken 0 0
[* , * , broken]:
maintain ignore :=
good 0 0
deteriorating 0 .5
broken .8 1 ;
param reward:
maintain ignore :=
good 1 2
deteriorating 1 2
broken -1 0 ;
param discount_factor := .8 ;
end ;
запускаемся: glpsol --math markov.mod -o output.txt :
Problem: markov
Rows: 7
Columns: 3
Non-zeros: 13
Status: OPTIMAL
Objective: total = 17.20538721 (MINimum)
No. Row name St Activity Lower bound Upper bound Marginal
------ ------------ -- ------------- ------------- ------------- -------------
1 total B 17.2054
2 bellman[good,maintain]
B 1.69697 1
3 bellman[good,ignore]
NL 2 2 7.76094
4 bellman[deteriorating,maintain]
NL 1 1 4.46128
5 bellman[deteriorating,ignore]
B 4.23906 2
6 bellman[broken,maintain]
NL -1 -1 2.77778
7 bellman[broken,ignore]
B 0.198653 -0
No. Column name St Activity Lower bound Upper bound Marginal
------ ------------ -- ------------- ------------- ------------- -------------
1 value[good] B 8.48485
2 value[deteriorating]
B 7.72727
3 value[broken]
B 0.993266
Karush-Kuhn-Tucker optimality conditions:
KKT.PE: max.abs.err = 0.00e+00 on row 0
max.rel.err = 0.00e+00 on row 0
High quality
KKT.PB: max.abs.err = 0.00e+00 on row 0
max.rel.err = 0.00e+00 on row 0
High quality
KKT.DE: max.abs.err = 8.88e-16 on column 1
max.rel.err = 2.22e-16 on column 3
High quality
KKT.DB: max.abs.err = 0.00e+00 on row 0
max.rel.err = 0.00e+00 on row 0
High quality
End of output
фабрика производит три вида продукции: P1, P2 и P3. известна цена на продукцию и спрос на каждый из видов продукции
| продукт | цена | спрос |
|---|---|---|
| P1 | 1200 | 35 |
| P2 | 2500 | 25 |
| P3 | 1400 | 30 |
на фабрике есть три цеха: W1, W2 и W3. для производства продукции P1 необходимы только технологические операции, производимые цехом W1, для P2 - W1 и W3, для производства P3 необходима обработка во всех трех цехах W1, W2 и W3. причем, если в цехах W1 и W2 продукция разных видов обрабатывается одинаково, и известна общая производительность этих цехов в единицах обработанной продукции, то в цехе W3 известны как временные затраты на обработку каждого из видов продукции, так и общий фонд рабочего времени сотрудников цеха
| W1, шт. | W2, шт. | W3, часов (P2/P3/фонд) |
|---|---|---|
| 40 | 20 | 8/2/80 |
из всех видов материалов, используемых при производстве продукции, ограниченным является только один, поставки его в неделю и потребности для каждого из видов продукции приведены в таблице
| поставки, кг | на ед. P1, кг | на ед. P2, кг | на ед. P3, кг |
|---|---|---|---|
| 50 | 0,8 | 0,6 | 0,7 |
необходимо составить производственный план, максимизирующий выручку от реализации продукции
файл модели, factory1.mod :
set P := 1..3 ;
set W := 1..3 ;
set W_Mass within W := {1, 2} ;
set W_Manual within W := W diff W_Mass ;
set Ops within P cross W := {(1, 1), (2, 1), (2, 3), (3, 1), (3, 2), (3, 3)} ;
set W_P {w in W} := setof {(i, w) in Ops} i ;
param price {p in P} ; # цена по видам продукции, руб.
param demand {p in P} ; # спрос в неделю по видам продукции, шт.
param cost {p in P} ; # расход материала на производство продукции, кг
param prod {w in W_Mass} ; # производительность цехов, шт.
param prod_m_cost {w in W_Manual, p in P} ; # затраты времени на единицу продукции, часы
param prod_m {w in W_Manual} ; # фонд рабочего времени, часы
param materials ; # поставки материалов в неделю, кг
var x {p in P} >= 0 ;
# максимальное количество изделий, которые могут пройти машинную обработку :
s.t.
mass {w in W_Mass} : sum {p in W_P[w]} x[p] <= prod[w] ;
# для цеха с ручной обработкой должно выполняться требование по рабочему времени :
s.t.
manual {w in W_Manual} : sum {p in W_P[w]} x[p] * prod_m_cost[w, p] <= prod_m[w] ;
# производство каждого из видов продукции не должно превосходить спрос :
s.t.
dem {p in P} : x[p] <= demand[p] ;
# затраты материала на производство не должны превышать поставок :
s.t.
mat : sum {p in P} x[p] <= materials ;
maximize
profit : sum {p in P} price[p] * x[p] ;
solve ;
end ;
# > glpsol -m factory1.mod -d factory1.dat -o factory1.sol файл данных, factory1.dat:
data ; # цена по видам продукции, руб. param price := 1 1200 2 2500 3 1400 ; # спрос по видам продукции, шт. param demand := 1 35 2 25 3 30 ; # расход материала на производство продукции, кг param cost := 1 0.8 2 0.6 3 0.7 ; # производительность цехов, шт. param prod := 1 40 2 20 ; # затраты времени на единицу продукции, часы param prod_m_cost default 0 := [3, 2] 8 [3, 3] 2 ; # фонд рабочего времени, часы param prod_m := 3 80 ; # поставки материалов, кг param materials := 50 ; end ;
запускаемся : glpsol -m factory1.mod -d factory1.dat -o factory1.sol
файл результатов, factory1.sol:
Problem: factory1
Rows: 8
Columns: 3
Non-zeros: 15
Status: OPTIMAL
Objective: profit = 61000 (MAXimum)
No. Row name St Activity Lower bound Upper bound Marginal
------ ------------ -- ------------- ------------- ------------- -------------
1 mass[1] NU 40 40 1200
2 mass[2] B 0 20
3 manual[3] NU 80 80 162.5
4 dem[1] B 30 35
5 dem[2] B 10 25
6 dem[3] B 0 30
7 mat B 40 50
8 profit B 61000
No. Column name St Activity Lower bound Upper bound Marginal
------ ------------ -- ------------- ------------- ------------- -------------
1 x[1] B 30 0
2 x[2] B 10 0
3 x[3] NL 0 0 -125
Karush-Kuhn-Tucker optimality conditions:
KKT.PE: max.abs.err = 0.00e+00 on row 0
max.rel.err = 0.00e+00 on row 0
High quality
KKT.PB: max.abs.err = 0.00e+00 on row 0
max.rel.err = 0.00e+00 on row 0
High quality
KKT.DE: max.abs.err = 0.00e+00 on column 0
max.rel.err = 0.00e+00 on column 0
High quality
KKT.DB: max.abs.err = 0.00e+00 on row 0
max.rel.err = 0.00e+00 on row 0
High quality
End of output оптимальное значение целевой функции (строка profit - по имени целевой функции, введенному в файле модели) 61000 достигается при производстве 30 штук продукции первого вида (значение из столбца Activity напротив переменной, соответствующей плану выпуска первой продукции - x[1]) и 10 штук продукции второго вида (аналогично, для переменной x[2]). Переменная x[3] помечена статусом NL (небазисная переменная, значение которой находится на нижней границе); действительно, в решении эта переменная принимает значение 0 (столбец Activity), а это является нижней границей области допустимых значений неотрицательной переменной. поскольку данная переменная является небазисной, она имеет ненулевую приведенную цену -125 (столбец Marginal). то есть, увеличение значения переменной x[3] на единицу повлечет за собой снижение значения целевой функции на 125 единиц
судя по значениям столбца St первой таблицы, из всех ограничений задачи активными являются два: mass[1] (ограничение производительности первого цеха) и manual[3] (ограничение производительности третьего цеха). статус NU (неравенство с активным ограничением сверху) согласуется с тем, что значение левой части ограничения (столбец Activity) для данных ограничений совпадает с верхней границей (столбец Upper bound). столбец Marginal содержит теневые цены для этих ограничений (1200 и 162,5 соответственно)
основные показатели, используемые при анализе чувствительности - теневую и приведенные цены - можно извлечь из таблицы результатов, но в GLPK предусмотрен специальный вид отчета по анализу чувствительности, содержащий расширенную информацию. в частности, в этом отчет можно использовать для определения интервалов осуществимости, то есть, границ, в которых могут изменяться правые части ограничений и коэффициенты целевой функции без изменения состава базисных переменных
для создания такого отчета при вызове glpsol следует указать параметр --ranges, задав имя файла, куда отчет должен быть помещен. например:
$> glpsol -m factory1.mod -d factory1.dat --ranges factory1_sens.sol -o factory1.sol
анализ чувствительности активных границ проводится только для активных ограничений. для каждого такого ограничения в ходе исследования активная граница (правая часть ограничения) изменяется как в сторону уменьшения, так и в сторону увеличения. каждому ограничению в отчете соответствует две строки: первая характеризует результат, полученный при уменьшении активной границы, вторая - при увеличении. изменение активной границы ведет к изменению значений переменных и продолжается до тех пор, пока какая-либо из базисных переменных не достигает границы области допустимых значений. соответствующее значение активной границы называется в GLPK break point и по сути является границей интервала осуществимости решения по отношению к исследуемому ограничению. в отчете по анализу чувствительности диапазон исследуемых значений границы приводится в столбце Activity range. соответствующие значения целевой функции указаны в столбце Obj value at break point, а символические имена ограничивающих базисных переменных - в столбце Limiting variable. если активную границу можно увеличивать или уменьшать бесконечно, то в столбце Activity range будет указано +Inf или –Inf соответственно, а столбец Limiting variable будет пуст
файл factory1_sens.sol:
GLPK 4.65 - SENSITIVITY ANALYSIS REPORT Page 1
Problem: factory1
Objective: profit = 61000 (MAXimum)
No. Row name St Activity Slack Lower bound Activity Obj coef Obj value at Limiting
Marginal Upper bound range range break point variable
---------- ------------- ------------- ------------- ------------- ------------- ------------
1 mass[1] NU 40.00000 . -Inf 10.00000 -1200.00000 25000.00000 x[1]
1200.00000 40.00000 45.00000 +Inf 67000.00000 dem[1]
2 mass[2] BS . 20.00000 -Inf . -Inf 61000.00000
. 20.00000 30.00000 125.00000 61000.00000 x[3]
3 manual[3] NU 80.00000 . -Inf 40.00000 -162.50000 54500.00000 dem[1]
162.50000 80.00000 200.00000 +Inf 80500.00000 dem[2]
4 dem[1] BS 30.00000 5.00000 -Inf 15.00000 -166.66667 56000.00000 x[3]
. 35.00000 40.00000 1300.00000 100000.00000 manual[3]
5 dem[2] BS 10.00000 15.00000 -Inf 5.00000 -500.00000 56000.00000 x[3]
. 25.00000 10.00000 +Inf +Inf
6 dem[3] BS . 30.00000 -Inf . -Inf 61000.00000
. 30.00000 20.00000 125.00000 61000.00000 x[3]
7 mat BS 40.00000 10.00000 -Inf 10.00000 -1200.00000 13000.00000 mass[1]
. 50.00000 40.00000 +Inf +Inf
8 profit BS 61000.00000 -61000.00000 -Inf 25000.00000 -1.00000 . mass[1]
. +Inf 61000.00000 +Inf +Inf
GLPK 4.65 - SENSITIVITY ANALYSIS REPORT Page 2
Problem: factory1
Objective: profit = 61000 (MAXimum)
No. Column name St Activity Obj coef Lower bound Activity Obj coef Obj value at Limiting
Marginal Upper bound range range break point variable
----- ------------- ------------- ------------- ------------- ------------- ------------
1 x[1] BS 30.00000 1200.00000 . 15.00000 1033.33333 56000.00000 x[3]
. +Inf 35.00000 2500.00000 100000.00000 manual[3]
2 x[2] BS 10.00000 2500.00000 . 5.00000 2000.00000 56000.00000 x[3]
. +Inf 10.00000 +Inf +Inf
3 x[3] NL . 1400.00000 . -6.66667 -Inf 61833.33333 dem[1]
-125.00000 +Inf 20.00000 1525.00000 58500.00000 mass[2]
End of reportнапример, ограничение mass[1] является активным в верхней границе (статус NU), по условию задачи равной 40 (Upper bound). соответственно, в ходе анализа чувствительности эта граница изменялась в сторону уменьшения до тех пор, пока при значении 10 переменная x[1] не приняла граничного значения (стала равной нулю), и в сторону увеличения до тех пор, пока при значении 45 не активизировалось ограничение dem[1] (спрос на продукцию первого вида). соответствующие значения целевой функции - 25000 и 67000
анализ чувствительности коэффициентов целевой функции при небазисных переменных является достаточно простым, поскольку изменение коэффициента целевой функции при небазисной переменной ведет к эквивалентному изменению ее приведенной цены. для каждой небазисной переменной коэффициент при этой переменной в целевой функции изменяется как в сторону уменьшения, так и в сторону увеличения. каждой переменной в отчете соответствует две строки: первая описывает результат, полученный при уменьшении коэффициента, вторая - при увеличении. изменение коэффициента ведет к изменению приведенной цены анализируемой переменной и влияет на приведенные цены других переменных. базисное решение остается допустимым и оптимальным до тех пор, пока приведенная стоимость сохраняет свой знак. в отчете по анализу чувствительности диапазон значений коэффициента выводится в столбце Obj coef range. если коэффициент можно увеличивать или уменьшать бесконечно без смены базиса, то в столбце Obj coef range будет указано +Inf или –Inf соответственно
например, переменная x[3] - небазисная переменная, принимающая в решении наименьшее допустимое значение (статус NL). коэффициент при этой переменной в исходной задаче равен 1400, а приведенная стоимость равна -125 (столбец Obj coef/Marginal). при уменьшении коэффициента при x[3] потенциальный вклад x[3] в целевую функцию становится еще меньше; что отражается и в значении -Inf в столбце Obj coef range. то есть, можно сколько угодно уменьшать данный коэффициент без изменения базиса
уменьшение цены на третий вид продукции не приведет к пересмотру производственного плана, в частности, к включению этого вида продукции в план. при увеличении же коэффициента до 1525 состав базисных переменных изменяется и разница между исходным значением коэффициента и пороговым равна 125, то есть, объяснимо совпадает с приведенной ценой переменной x[3]
анализ чувствительности коэффициентов целевой функции при базисных переменных - для каждой базисной переменной коэффициент при этой переменной в целевой функции изменяется как в сторону уменьшения, так и в сторону увеличения. каждой переменной в отчете соответствует две строки: первая описывает результат, полученный при уменьшении коэффициента, вторая - при увеличении. изменение коэффициента ведет к изменению приведенной цены небазисных переменных. при определенном изменении какая-то из цен становится равной нулю, а это означает, что даже малейшее дальнейшее изменение приведет к включению в базис соответствующей переменной. в отчете по анализу чувствительности диапазон значений коэффициента выводится в столбце Obj coef range. соответствующие значения целевой функции приведены в столбце Obj value at break point, символические имена ограничивающих небазисных переменных - в столбце Limiting variable, а значения исследуемой переменной, которые она принимает на концах диапазона - в столбце Activity range. если коэффициент можно увеличивать или уменьшать бесконечно без смены базиса, то в столбце Obj coef range будет указано +Inf или –Inf соответственно
например, переменная x[1] является базисной и принимает в решении значение 30 (столбец Activity), а коэффициент в целевой функции при этой переменной равен 1200 (столбец Obj coef/Marginal). анализ показал, что переменная остается базисной при коэффициенте в целевой функции не ниже 1033.33333, при дальнейшем снижении коэффициента x[1] выйдет из базиса и ее место займет x[3]. значение переменной x[1] при граничном значении коэффициента будет равно 15 (столбец Activity range), а значение целевой функции – 56000 (столбец Obj value at break point). увеличивать коэффициент при x[1] без смены базиса можно вплоть до значения 2500. значения переменных при этом будут изменяться, значение целевой функции возрастет до 100000, но состав ненулевых переменных и активных ограничений будет неизменен. дальнейшее повышение коэффициента ограничено manual[3]
ограничение manual[3] в решении является активным. это ограничение является ограничением нестрогого неравенства, а значит, при переходе к канонической форме в него должна быть добавлена дополнительная переменная, значение которой равно разнице между левой и правой частями исходного ограничения. при активном ограничении левая и правая части равны, а значит, дополнительная переменная не содержится в базисе. при повышении коэффициента целевой функции при x[1] выше 2500 дополнительная переменная входит в базис, то есть, левая его часть соответствующего ограничения становится строго меньше правой. это означает, что при ценах на продукцию первого вида выше 2500 план производства изменяется таким образом (в пользу продукции первого вида), что третий цех становится недогруженным