i-808190579

Если исходная точка x X , то	x	x rnd(1) x x . Процедура
0	0i	i	i	i
повторяется N1 раз. Это обеспечивает вписывание исходной точки х0 в огра-
ничения задачи.
Во всех вершинах исходного симплекса найти значение
Qj Q x j	, j ,...,n ,

и если текущая j-я вершина симплекса не удовлетворяет ограничениям задачи (8.1) Qj 900000 j при поиске максимума и Qj 900000 j при поис-

ке минимума.

Положить номер шага поиска k n ; Yk Qk – массив для построения графика Qk ; Х – массив координат вершин симплекса;

2) из всех вершин текущего симплекса, выбрать вершину с наименьшим (поиск максимума) с наибольшим (поиск минимума) значением Q j ,

причем j p, где р – номер новой вершины на предыдущем шаге поиска (запрет возврата)

Qs min Q j – поиск максимум;

j p

Qs max Q j – поиск минимума.

j p

Положить s j , p s , где s – номер отраженной и новой вершины симплекса;

3) вычислить координаты вершины Vsн нового симплекса по формуле

xsiн	n	n	xsi ,i ,...,n ;
	x ji	n
	n j

4) проверить выполнение ограничений (8.1) для точки с координатами

xsiн :

если эта точка удовлетворяет ограничениям, то перейти к п.5,

		иначе Qн 90000 k и	перейти к п.6 (при		поиске	минимума
		s
Qн 90000 k );
s
	5)	определить значение Qн	в вершине с координатами			xн , i ,n ;
		s				si
Y Qн		– массив для построения графика Q		;
k	s			k
	6) положить k k , m j m j ; ms			, здесь m j	– число последова-

тельных шагов поиска, в которых вершина с номером j не отражалась. Для новой (отраженной) вершины ms ;

7) присвоить	j s , x	ji	xн	,i ,...,n ,	Q	j	Qн . Вычислить координаты
			si				s

остальных вершин текущего симплекса

	xн	xн	x		xн		k	, j ,...,n , j s, i ,...,n ,
				ji
	ji	si			si
							k

131

где γ	k	e μk ;
	8) провести анализ m j , j ,...,n :
			если m j n , перейти к п.9);
			иначе d j , вычислить	Q Q xн (повторный эксперимент),
				d	d

md , перейти к п.9);

9) проверить число шагов поиска k:

если k N , перейти к п.2);

если k N , перейти к п.10);

10) вывести на печать оптимальные значения x*j,i , Q*j , i ,...,n, j ,...,n и траекторию движения координат вершин симплекса Х k для двух выбранных переменных xr , xe на двухкоординатный график и график Yk от k.

Для cистемы Matlab.

Программа состоит из следующих файлов:

lab8.m;

Q.m;

fogr.m;

SimplexSearch.m;

CalculateEdges.m;

CheckPosConstraints.m;

SimplexSearchBody.m.

Вданном списке полужирным шрифтом выделены файлы, в которых производится настройка поиска.

Файл lab8.m

Вфайле lab8.m происходит инициализация поиска, вывод выходных данных и построение графиков истории поиска. Здесь необходимо задать следующие параметры:

kp – режим поиска – поиск минимума (значение 0) или максимума (значение 1);

n – размерность пространства поиска;

X0 – начальное положение поиска (массив переменных x1, x2, …, xn);

xmin, xmax – позиционные ограничения (массивы);

axes – номера переменных, для которых нужно построить графики истории поиска (массив из двух чисел);

gridRes – разрешение сетки графиков истории поиска.

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

132

Файл Q.m

В файле Q.m задается целевая функция. Переменные передаются в массиве x; для обращения к конкретной переменной необходимо использовать синтаксис x(i), где i – порядковый номер переменной. Результат записывается в переменную ret.

Файл fogr.m

Вфайле fogr.m описываются функциональные ограничения. Аналогично с Q.m, переменные передаются в массиве x. Если данное значение x находится вне зоны поиска, то в ret записывается положительное число, в противном случае в ret записывается отрицательное число или ноль.

Файл SimplexSearch.m

Файл SimplexSearch.m является основным файлом программы поиска: функция SimplexSearch принимает параметры поиска – режим поиска, количество переменных, целевую функцию и др., и возвращает историю поиска – массив всех точек траектории поиска и значений целевой функции в них.

Вданном файле находятся следующие параметры тонкой настройки

поиска:

nu – точность поиска;

N1 – число итераций для «вписывания» исходной точки Х0 симплекса в ограничения задачи.

Файл CalculateEdges.m

Функция CalculateEdges возвращает длины ребер начального и конечного симплексов, а также расстояние между точками xmin и xmax (см lab8.m).

Файл CheckPosConstraints.m

В этом файле находится функция, принимающая массив переменных поиска x и возвращающая количество тех из них, для которых были нарушены позиционные ограничения.

Файл SimplexSearchBody.m

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

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

1.В соответствии с вариантом W выберите модель, функциональное ограничение и исходную точку поиска из лабораторной работы № 7. Позиционные ограничения сформируйте самостоятельно, исходя из графического представления целевой функции, приведенного в лабораторной работе № 7.

2.Подпишите карандашом комментарии на листинге программы, согласно алгоритма.

3.Внесите изменения в программу (Mathcad): n:= , kp:= , xmin:= ,

xmax:= , x0:= , a:= , b:= , xx:= если требуется по модели, Q(y):= , q1(y):= , r:= , e:= .

133

4. Постройте также графики функций Q(xr , xe ) и q1(xr , xe ) . Выделите карандашом допустимое множество Х.

5.Решите задачу стандартной процедурой Mathcad. Сопоставьте полученные результаты.

6.Используя график функции q1 (x) измените позиционное ограниче-

ние так, чтобы изменилась точка оптимума x* и повторите процесс поиска.

7.Решите задачу в системе Matlab, внеся соответствующие изменения

впрограмме.

Выводы по работе должны включать:

1)сравнительный анализ полученных результатов с градиентным методом и методом системы Mathcad, Matlab;

2)объяснение основных переменных программы.

КОНТРОЛЬНЫЕ ВОПРОСЫ

1. Что представляет собой симплекс для n 1 и

n 3 ?

2.Зачем уменьшают размер симплекса в процессе поиска?

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

4.Почему ПСМ называют методом с совмещенными пробными и рабочими шагами?

5.Какие достоинства и недостатки присущи ПСМ.

ЛИТЕРАТУРА

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

2.Дамбраускас А.П. Симплексный поиск. – М.: Энергия, 1979.-175 с.

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

4.Макаров Е. Инженерные расчеты в Mathcad 15: Учебный курс. –

СПб.: Питер, 2011. – 400 с.: ил.

134

Тестовый пример в системе Mathcad

Пос ледовательный с имплекс ный метод с пере менным шагом, запретом возврата, рас четом параметров поис ка. Поис к макс имума kp=1, минимума kp=0

n 2	- размернос ть вектора х 0.01 - точнос ть отыс кания экс тремума
N1 5	чис ло итерации для "впис ывания" ис ходной то чкив ограничениях задачи
	0
Позиционные ограничения OR IGIN 1 kp 0

	10						0
xmax		xmin
	10						0																						a0 2
																		2.8
																		2.0
	8																											1
														a					0.5					xx
X0	- ис ходная точка																												- параметры модели
	8																			1								1
																			0.5
1								параметры функциональных ограничений задачи
b
1
Q(y) a0		a	1	y		1		a		2	y		2	a	3	y			1	y		2	a	4	y			2	a	5	y 2		целевая функция
																												1				2
	n		bi yi
g1(y)											b1 4							- функциональные ограничения
	i 1						r 1							e 2					- координатные ос и для пос троения графика
Q1(x1 x2) a0			a		1		x1				a	2	x2 a				3		x1 x2 a							4	(x1)2 a				5	(x2)2

q(x1 x2) b1 x1 b2 x2 b1

Q1	q

135

Рас чет ис ходного и конечного ребра с имплекс а

( L0 rmax)

for

i 1 n

d

i

xmax xmin

i

i

rmin min(d)

1

L0

rmin

n (n 1) 2

2 n

2

1

n

di

2

rmax

2

i 1

( L0

rmax)

				1
LN		n		2
	2
		n 1

LN 0.012

L0 4.33 rmax 14.142

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

poz(y)

sum 0

for

i 1 n

sp

i

0

if y

i

xmax y

i

xmin

i

i

sp i 1

otherwise

sum sum sp i

sm sum

sm

Ос новной блок программы

( X Y )

k n 1

N

8

rmax

L0

1

(ln(L0) ln(LN))

N

(k) exp( k)

p n 2

w n 1

for t 1 N1

for

j 1 n 1

for i 1 n

0

if i j 1

1

if i j 1

2 i

(i

1)

i

if i

j 1

2(i

i

1)

i

i

i

X0 xmin rnd(1)

xmax xmin

xj i X0i L0

Xj i xj i

136

zzi

xj i

XT XT

fogr (y)		sp1 0	if g1(y) b1
		sp1 1	otherwise
		sp1

if poz(X0) 1	fogr (X0) 1

i

if i

j 1

1)

i

2(i

i

i

i

X0

xmin rnd(1)

xmax xmin

if poz(X0) 1 fogr (X0)

1

xj i X0i L0

Xj i xj i

zzi xj i

XT XT

for

j 1 n 1

z j XT j

f

j

j

Q z

j

j

f

j

900000 k

if

poz z

1 fogr

z

1

kp

1

if poz z j 1 fogr z j 1

f

900000 k

kp

0

j

Yj fj

mj 0

while k N

ft sort (f)

ft reverse (ft )

if kp 0

fm ft

1

for j 1 n 1

s j

if fj

fm

fm ft

2

if s

p

for j 1 n 1

s j

if fj

fm

p s

for i 1 n

n 1

2

(2 n)

xnovs i

xj i

xs i

n

n

j 1

znovi xnovs i

fs Q(znov)

fs 900000 k

if

(poz(znov) 1 fogr (znov) 1) kp

1

fs 900000 k

if (poz(znov) 1 fogr (znov) 1) kp

0

for j 1 n 1

mj mj

1

ms 0

k k 1

for i 1 n

xs i xnovs i

for j 1 n 1

i 1 n

for

x

xnov

x

xnov

137

(k)

if j s

(k 1)

j

i

s

i

j i

s i

for j 1 n 1

for

i 1 n

					ms 0
					k k 1
					for		i 1 n
					xs i xnovs i
					for		j 1 n 1
								i 1 n
					for
										xnov						x	xnov		(k)		if j s
						x
																			(k
							j	i				s		i		j i	s i			1)
							j 1 n 1
					for
					for			i 1 n
						Xw j i xj i
					Yk fs
					w w n 1
					mp max(m)
					( continue ) if							mp n 2
					for		j 1 n 1
					d j if mj									mp
					for		i 1 n
					zz x
						i				d i
					fd Q(zz)
					k k 1
					Y f
										d
					k
					md 0
				( X Y )
																	Результаты поис ка
																			1
			1							2
																		1	48.683
	1		5.835							6.75
																		2	9·105
	2		10.165							6.75
																		3	9·105
	3		8							10.5
																		4	54.1
	4		6.274						5.989
																		5	15.277
	5		9.726						5.989
																		6	15.444
	6		8							3
																		7	3.216
	7		5.924						5.382								Y
X																		8	3.352
	8		4.549							3
																		9	9·105
	9		7.3							3
																		10	9·105
	10		5.366						5.382
																		11	2.155
	11		4.27						3.483
																		12	1.564
	12		3.173						5.382
																		13	4.021
	13		2.077						3.483
																		14	3.741
	14		3.825						3.483
																		15	2.147
	15		2.951						4.997
																		16		...
	16		2.254							...

138

m 1 rows(X)		k 1 rows(Y)
Решение процедурой Mathcad
	y X0
	Given		y xmin		y xmax
		g1(y) b1
	xy MinimizeQ( y)			1.4		Q(xy) 0.04
			xy
				2.6
			z 0 10
	Qo Q(xy)			b1 b1 z
			x2(z)		b2
Перемещение с имплекс а в прос транс тве коор динат
10
8
Xm e 6
xy2
x2(z) 4
2
0	0	2	4	6	8	10
			Xm r z

Изменение целевой функции в процес с е поис ка

40

20

Yk

Qo

0

20

10

20

30

k

139

Тестовый пример в системе Matlab

%Режим поиска: 0 - минимум, 1 - максимум (задать) kp = 0;

%Размерность пространства (задать) n = 2;

%Начальная точка поиска (задать)

X0 = [8 8];

%Позиционные ограничения (задать) xmin = [0 0];

xmax = [10 10];

%Номера координатных осей для построения графиков (задать) axes = [1 2];

%Разрешение сетки графиков (задать) gridRes = 20;

%Вызов процедуры поиска

[X, Y] = SimplexSearch(kp, n, @Q, X0, xmin, xmax, @fogr);

%Вывод данных

%disp('Траектория поиска:');

%disp (X);

disp('Конечный симплекс:'); for i = 1 : n + 1

disp(['Координаты: ', num2str(X(end - i, :), '%f ')]);

disp(['Значение целевой функции: ', num2str(Q(X(end - i, :)), '%f ')]);

end

%Сетка для построения графиков

[x1, x2] = meshgrid(xmin(axes(1)) : ...

(xmax(axes(1)) - xmin(axes(1))) / gridRes : ...

xmax(axes(1)), ...

xmin(axes(2)) : ...

(xmax(axes(2)) - xmin(axes(2))) / gridRes : ...

xmax(axes(2)));

%Выбор первой части окна для построения графиков subplot(1, 2, 1)

axis([xmin(axes(1)) xmax(axes(1)) xmin(axes(2)) xmax(axes(2))]);

%Построение линий уровня целевой функции(двухмерный вид) z = zeros(length(x1), length(x2));

for i = 1 : length(x1)

for j = 1 : length(x2)

z(i, j) = Q([x1(i, j) x2(i, j)]);

end

end

contour(x1, x2, z, 'ShowText', 'on'); hold on;

%Построение нулевой линии уровня функции ограничения (двухмерный вид) z = zeros(length(x1), length(x2));

for i = 1 : length(x1)

for j = 1 : length(x2)

z(i, j) = fogr([x1(i, j) x2(i, j)]);

end

end

contour(x1, x2, z, [0 0], 'LineWidth', 2);

%Построение траектории поиска (двухмерный вид) plot(X(:, axes(1)), X(:, axes(2)), '.-black'); hold off;

%Выбор второй части окна для построения графиков

140