i-808190579

2.На что влияет коэффициент α ?

3.Покажите графически как учитывают позиционные ограничения.

4.В чем заключается случайная составляющая комплекс-метода?

ЛИТЕРАТУРА

1.Масальский, Г.Б. Математические основы кибернетики: учебное пособие / Г.Б. Масальский. – Красноярск: Сиб. федер. ун-т, 2014. – 215 с. (электронная версия).

2.Банди, Б. Методы оптимизации. Вводный курс: Пер. с англ. / Б. Банди. – М.: Радио и связь, 1988. – 128 с., ил.

5.Пантелеев, А. В. Методы оптимизации в примерах и задачах: Учеб. пособие /А.В. Пантелеев, Т.А. Летова. – М.: Высш. шк., 2002. – 544 с.: ил.

6.Кетков Ю.Л., Кетков А.Ю., Шульц М.М. MATLAB 7: программирование, численные методы. – Спб.: БВХ-Петербург, 2005. – 752с.: ил.

151

Листинг программы

%*************************Комплекс метод***********************************

%Примечание: позиционные, функциональные ограничения задаются в

%соответствующих процедурах (); условие останова выбирается оператором

%вручную в тексте программы перед циклом поиска.

%%Очистка переменных clc; clear;

%%Определение параметров

%------------------------				Основные параметры--------------------------------
n =	2 ;			% размерность поиска
L =	2*n	;		% число вершин комплекса
alf	= 1.3 ;			% коэффициент растяжения комплекса
kp = 0;				%поиск минимума
%kp	= 1		%поиск максимума
%Позиционные ограничения
	xmin = [0			0];
	xmax = [10 10];
%-------------------------				Условие останова---------------------------------
N =	100	;	%	количество шагов поиска

SKOx = 0.01; SKOq = 0.001;

%СКО отклонения координат вершин текущего комплекса и значений

%целевой функции в них соответственно

%------------

Параметры для вывода графических данных на экран--------------

%Одновременно на график можно вывести результат поиска для двух

%переменных, если n > 2, то необходимо определить какие это будут

%переменные

r = 1	;	% абсцисса для	графика
e = 2	;	% ордината для	графика
%--------------		Служебные	параметры		(не редактируются)----------------------
k = 1;		% счетчик числа	шагов	поиска
x = 0;		% определение массива		для	хранения текущих вершин комплекса
q = 0;		% определение вектора		для	хранения значения функции в текущих
		% вершинах комплекса
X = 0;		% определение массива		для	хранения всей истории вершин комплекса
Q = 0;		% определение вектора		для	хранения всей истории значений функции

%в вершинах комплекса

%%Формирование исходного комплекса

%формирование начальной (исходной) точки X1(X0),удовлетворяющей %ограничениям

for i = 1 : n

x0(1,i) =xmin(i)+rand(1)*(xmax(i)-xmin(i));

end

while ogr(x0)>=0 % задать по условию задачи >= или <= for i=1:n

x0(1,i)=x0(1,i)+ (-1+rand(1))*(xmax(i)-xmin(i));

end

end

for i=1:n x(1,i)=x0(1,i); X(1,i)=x(1,i);

end

q(1)=funct(x,n,1);

Q(k)=q(1);

%формирование остальных вершин комплекса for j = 2 : L

k = k + 1; % счетчик итераций в исходном комплексе for i = 1 : n

% построение текущей j-ой точки исходного комплекса

152

x(j,i) = xmin(i)+ (-1+rand(1))*(xmax(i)-xmin(i)); xr(1,i)=x(j,i);

end

xr=prov_p(xr,n,xmin,xmax); xr=prov_q(xr,x0,n);

%формирование центра тяжести j вершин,удовлетворяющих ограничениям for m=1:n

x0(1,m) = sum(x(:,m))/j ;

end

for i=1:n

% запись в архив (для графика)

X(k,i) = x(j,i);

end

%вычисление значения целевой функции в текущей точке q(j) = funct(x,n,j);

%запись в архив (для графика

Q(k) = q(j);

end k1=k;

Sx = SKO_x(x,L,n); % вычисление СКО точек текущего комплекса SQ = SKO_Q(q,L); % вычисление СКО функции текущего комплекса %% Поиск

%Для останова алгоритма поиска необходимо выбрать одно из условий:

%1) определенное количество шагов;

%2) факт сжатия вершин в достаточно малый объем, определяемый через СКО

%точек и значений функции вершин текущего комплекса.

%Для этого необходимо снять комментарий с нужной строки, и наоборот -

%закомментировать другую (может быть выбрано только одно условие !!!).

% for	k = k1 :N	%	закомментировать	по условию 2
while	SQ > SKOq	&& Sx >	SKOx k = k + 1;	% закомментировать по условию 1

%вызов процедуры определения индекса худшей вершины s = sortirovka(q,L,kp);

%запомним значение отклика в худшей вершины комплекса q_st = q(s);

%Отражение худшей точки комплекса

for i = 1 : n

%центр тяжести лучших (L-1) вершин комплекса x0(1,i) = (sum(x(:,i))-x(s,i))/(L-1);

%новая отраженная точка комплекса

xr(1,i) = (1+alf)*x0(1,i)-alf*x(s,i);

end

%Для того, чтобы исключить ограничения из задачи поиска, необходимо

%закомментировать нужные строки, расположенные ниже % вызов процедуры проверки позиционных ограничений;

%при отсутствии закомментировать

xr = prov_p(xr,n,xmin,xmax);

%вызов процедуры проверки функциональных ограничений;

%при отсутствии закомментировать

xr = prov_q(xr,x0,n);

q_nov = funct(xr,n,1); % значение целевой функции в новой точке

%Проверка шага на успех

%если значение целевой функции в новой точке хуже, чем в старой, то

%для поиска минимума

if kp	== 0
while	q_nov>=q_st
for i	= 1 : n

%сдвигаем новую точку на половину расстояния между ней

%и центром тяжести лучших вершин комплекса

xr(1,i) = (xr(1,i)+x0(1,i))/2;

153

end xr=prov_p(xr,n,xmin,xmax); xr=prov_q(xr,x0,n);

%значение целевой функции в пересчитанной точке q_nov = funct(xr,n,1);

end

%для поиска максимума

elseif kp == 1 while q_nov<=q_st

for i = 1 : n

%сдвигаем новую точку на половину расстояния между ней

%и центром тяжести лучших вершин комплекса

xr(1,i) = (xr(1,i)+x0(1,i))/2; end

xr=prov_p(xr,n,xmin,xmax); xr=prov_q(xr,x0,n);

% значение целевой функции в пересчитанной точке q_nov = funct(xr,n,1);

end end

% Сохранение проверенной вершины for i = 1 : n

x(s,i) = xr(1,i); % новая точка комплекса X(k,i) = x(s,i); % запись в архив (для графика)

end

q(s) = q_nov; % запишем в комплекс значение целевой функции % в новой точке

Q(k) = q(s); % запись в архив (для графика)

Sx = SKO_x(x,L,n); % вычисление СКО точек текущего комплекса SQ = SKO_Q(q,L); % вычисление СКО функции текущего комплекса

end

%%Построение графика grafic(X,x,Q,q,r,e,xmin,xmax)

%%Вывод результатов поиска disp('Последний комплекс') disp(x)

disp('Значения функции в вершинах последнего комплекса') disp(q)

disp('Количество шагов поиска') disp(k)

function Q = funct(x,n,j)

a0 = 2; a = [-2.8 -2.0 0.5 1 0.5];

Q=a0+a(1)*x(j,1)+a(2)*x(j,2)+a(3)*x(j,1)*x(j,2)+a(4)*x(j,1)^2+a(5)*x(j,2)^2;

return

function grafic(X,x,Q,q,r,e,xmin,xmax) %-------------------------------------------------------------------------

% Процедура построения графиков

%-------------------------------------------------------------------------

[XX,YY] = meshgrid(xmin(r):xmax(r), xmin(e):xmax(e));

154

%задаем сетку значений переменных, желательно равностороннюю

[s1, s2] = size(XX);

%определяем размерность

for g = 1 : s2

for j = 1 : s1

xc = [XX(j,g) YY(j,g)];

%вычисляем значение целевой функции во всех узлах сетки

Z(j,g) = funct(xc,2,1);

%вычисляем значение функции ограничения во всех узлах сетки

Z1(j,g) = ogr(xc);

end

end

%-------------------------------------------------------------------------

% Двухмерный график с линиями равного уровня

%-------------------------------------------------------------------------

figure(1) hold on grid on

plot(X(:,r),X(:,e),'-*m')

% вся история поиска plot(x(:,r),x(:,e),'-*k')

%последний комплекс (окрашен в другой цвет) contour(XX,YY,Z,30)

%линии равного уровня

hold off xlabel('x1') ylabel('x2') axis equal

%-------------------------------------------------------------------------

% Трехмерный график

%-------------------------------------------------------------------------

% figure(2) % hold on

% surf(XX,YY,Z);

% % построение модели целевой функции

% surf(XX,YY,Z1);

% % построение модели ограничения

% colormap jet

% % цветовая палитра

% shading interp;

% % сглаживание цвета

% colorbar

% % шкала цветовых оттенков

% grid on

% plot3(X(:,r),X(:,e),Q(:),'-*k') % plot3(x(:,r),x(:,e),q(:),'-*w')

% % нанесение истории поиска с последним комплексом

% hold off

function q_ogr = ogr(x)

% Функциональное ограничение b1=4;

q_ogr = -1*x(1,1)-1*x(1,2)+ b1; % программируется оператором return

155

function x_nov = prov_p(xr,n,xmin,xmax)

%Процедура проверки позиционных ограничений

%Определение минимальных и максимальных ограничений для каждой из

%переменных. Если необходимо, можно задавать ограничения отдельно для

%каждой из координат.

delta = 0.001; % величина приращения при выходе за границу

for i = 1:n

while xr(1,i) <= xmin(i) xr(1,i) = xr(1,i) + delta;

end

while xr(1,i) >= xmax(i) xr(1,i) = xr(1,i) - delta;

end

end

x_nov = xr;

return

function x_nov = prov_q(x,x0,n)

%Процедура проверки функциональных ограничений

%вычисление функции ограничения в текущей точке

% условие задается оператором в зависимости от ограничения m=0;

while ogr(x) >= 0 && m<=100 for i = 1 : n

x(1,i) = (x(1,i)+x0(1,i))/2;

%сдвигаем точку на половину расстояния до центра тяжести лучших

%вершин комплекса

end m=m+1;

end

x_nov = x; return

function SQ = SKO_Q(q,L)

% Процедура вычисления СКО значений целевой функции в вершинах комплекса

Qsr = sum(q(:))/L; Qc = 0;

for j = 1 : L

Qs(j) = (q(j)-Qsr)^2; Qc = Qc+Qs(j);

end

SQ = sqrt(Qc/(L-1));

return

function Sx = SKO_x(x,L,n)

% Процедура вычисления СКО вершин комплекса for i = 1:n

156

xc = 0;

xsr(i) = sum(x(:,i))/L;

for j=1:L

xc = xc+(x(j,i)-xsr(i))^2;

end

Sx2(i) = xc/(L-1);

end

Sx = sqrt(sum(Sx2(:))/n);

return

function s = sortirovka(q,L,kp)

% Процедура определения индекса худшей точки. s = 1;

if kp == 0

maxQ = q(1); for j = 2 : L

if q(j) > maxQ maxQ = q(j); s = j;

end

end

elseif kp == 1

minQ = q(1);

for j = 2 : L

if q(j) < minQ minQ = q(j); s = j;

end

end

end

return

Последний комплекс

1.3709 2.6301

1.3847 2.6172

1.3881 2.6126

1.3952 2.6052

157

Значения функции в вершинах последнего комплекса

0.0421 0.0427 0.0409 0.0405

Количество шагов поиска

49

158

ЛАБОРАТОРНАЯ РАБОТА № 10

ЛИНЕЙНОЕ ПРОГРАММИРОВАНИЕ

Цель работы – разработка производственного плана с использованием методов линейного программирования (ЛП).

КРАТКОЕ ТЕОРЕТИЧЕСКОЕ ОПИСАНИЕ

Теоретические сведения приведены в разделе 4 [1], а также в [2, 3, 4]. Постановка задачи линейного программирования в матричной форме

для общего случая имеет вид:

Q(x) cx max ,

x X

где X x : x En ,x 0, Ax b .

Исходной задаче линейного программирования соответствует двойственная задача вида:

G(λ) bT λ min ,

λ

где λ : λ E ,λ 0,A λ c – вектор-столбец двойственных переменных.

m T

Для решения задач линейного программирования в Matlab использовать функцию linprog.

Обращение к функции

[x, Q]=linprog (с, A, b, A1, b1, xmin, xmax).

ЗАДАНИЕ И ПОРЯДОК ВЫПОЛНЕНИЯ РАБОТЫ

Требуется составить оптимальный по реализации производственный план выпуска двух видов изделий в количестве x1 и x2 при определенных

возможностях четырех видов машин. Оба изделия последовательно обрабатываются этими машинами в течение aij времени на каждом i-ом станке. При

этом ежедневно каждый вид машин может работать не более bi часов.

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

В соответствии с условиями задачи (таблица 1) и вариантом (таблица 2) необходимо:

1)составить задачу линейного программирования, т.е. записать вид целевой функции и ограничений;

2)составить двойственную ей задачу линейного программирования;

3)решить графически (вручную) исходную задачу линейного программирования;

159

4) решить исходную и двойственную задачу ЛП средствами Mathcad и Matlab;

Таблица 1

Условия задачи линейного программирования

Вид машин							Вид изделия								Допустимое время
					x1						x2				работы машин
1-й					a						a					b1
					11						12
2-й					a	21					a	22				b2
3-й					a	31					a	32				b3
4-й					a	41					a	42				b4
Доход от реали-					с1						с2
зации изделия
																Таблица 2
Исходные данные для решения задачи линейного программирования
										Вариант
Коэффициент
		1		2		3		4		5		6	7		8	9	10
a		1		6		2		1		3		1	2		0	3	0,5
11
a		1		0		0		1		1		1	0		1	2	1,5
12
b1		18		36		40		18		18		18	10		18	24	12
a		0,5		1		1		0,5		1		1	0		2	0	1
21
a		1		2		1		1		1		0	2		1	2	0,5
22
b2		12		24		24		12		10		10	12		24	16	12
a		1		2		0		1		3		0	3		1	1	0
31
a		0		0		1		0		0		1	1		0	0	2
32
b3		14		24		16		16		15		14	24		8	7	12
a		0		0		1		0		0		1	1		0,5	2	1
41
a		1		2		2		1		1		0	1		0,5	2	0
42
b4		9		18		33		9		8		8	10		10	18	11
с1		4		4		1		6		3		1	4		1	2	1
с2		6		6		2		4		2		2	6		5	3	1
							160