Федеральное государственное бюджетное образовательное учреждение
высшего профессионального образования
«Ижевский государственный технический университет»
им. М.Т.Калашникова
Факультет «Информатика и вычислительная техника»
Кафедра «Вычислительная техника»
Пояснительная записка
к курсовой работе по дисциплине «Методы оптимизации»
на тему «Численные методы оптимизации функции нескольких переменных»
Вариант 12
Выполнил: студент гр.8-78-2
Принял: к.т.н., доцент
Ижевск 2014
Оглавление
1. Введение 3
2. Постановка задачи 4
3. Решение задачи с помощью аналитических методов 5
4. Метод Хука-Дживса с использованием метода Пауэлла 7
4.1. Метод Хука-Дживса 7
4.2. Метод Пауэлла 9
4.3. Описание программы 11
4.4. Результат выполнения программы 13
5. Метод наискорейшего спуска с использованием метода Фибоначчи 17
5.1. Метод наискорейшего спуска 17
5.2. Метод Фибоначчи 19
5.3. Определение отрезка унимодальности 21
5.4. Описание программы 22
5.5. Описание работы программы 24
6. Выводы 25
7. Приложение. Программа на ЯВУ 26
8. Список используемых источников 30
Введение
Основной целью выполнения настоящей курсовой работы является применение на практике основных численных методов оптимизации функций нескольких переменных.
Оптимизация функции представляет собой поиск минимального значения функции на некотором интервале значений аргументов функции.
Существует ряд методов оптимизации.
Их можно разделить на прямые методы, методы первого порядка, методы второго порядка [1].
Прямые методы поиска или методы нулевого порядка это методы, в которых при поиске экстремума используется информация только о самой функции и не используется информация о ее производных. Плюсом таких методов является возможность оптимизации функций, аналитическое представление которых неизвестно, т.е. эти функции определены только алгоритмически. К прямым методам оптимизации относятся метод Гаусса, метод Хука и Дживса, метод Пауэлла и другие.
Методы первого порядка при поиске решения используют не только информацию о самой функции, но и о ее производных первого порядка. К таким методам относятся различные градиентные алгоритмы, такие как метод наискорейшего спуска, метод сопряженных градиентов и т.д.
Методы второго порядка при поиске решения используют информацию о самой функции и о еѐ производных первого и второго порядка. Сюда относятся метод Ньютона и его модификации.
Некоторые из этих методов будем использовать в данной курсовой работе в применении к функции двух переменных f(x,y).
Постановка задачи
Построить график функции f(x,y) = 3x2+2y2+6x+3y, используя любой математический пакет.
Решить задачу безусловной оптимизации следующими способами:
Найти стационарные точки аналитическим методом, исследовать функцию в этих точках.
Найти точки экстремумов функции методом Хука-Дживса с использованием метода Пауэлла.
Найти точки экстремумов функции градиентным методом наискорейшего спуска с использованием метода Фибоначчи.
Решение задачи с помощью аналитических методов
Исходная функция f(x,y) = 3x2+2y2+6x+3y всюду непрерывна и дифференцируема.
Найдем стационарные точки функции. Для этого найдем градиент этой функции:
∇f(x,y)=(■(6x+6@4y+3))
и приравняем ее составляющие нулю:
{■(6x+6=0,@4y+3=0.)=>{■(x=-1,@y=-0.75.)┤ ┤
Исходная функция является суммой выпуклых функций, поэтому и она сама является выпуклой функцией. Следовательно, найденная стационарная функция является точкой минимума функции.
Найдем значение функции в точке минимума:
f(-1,-0.75) = 3*(-1)2+2*(-0.75)2+6(-1)+3(-0.75) = -4.125
С помощью математического пакета MathCad построим график исходной функции, тем самым визуально проверим, правильно ли была найдена точка минимума (смотрите рисунок 1).
Рисунок 1 – График функции f(x,y) = 3x2+2y2+6x+3y, построенный с помощью математического пакета MathCad
Как видно из рисунка 1, аналитическое нахождение точки минимума функции f(x,y) = 3x2+2y2+6x+3y дало правильный результат.
Кроме того, отметим, что функция имеет ровно одну точку минимума.
Метод Хука-Дживса с использованием метода Пауэлла
Метод Хука-Дживса
Другое название функции – метод конфигураций.
В данном алгоритме предлагается логически простая стратегия поиска, в которой используются априорные сведения о топологии функции и, в то же время, отвергается уже устаревшая информация об этой функции. Алгоритм включает два основных этапа [1]:
1)исследующий поиск вокруг базисной точки x ̅^k;
2)поиск по "образцу", т.е. в направлении, выбранном для минимизации.
В первую очередь задается начальная точка поиска x ̅^0 и начальное приращение (шаг) 〖∆x ̅〗^0. После этого начинается исследующий поиск.
Исследующий поиск. Делаем пробный шаг по переменной x_1, т.е. определяем точку x_1^0+∆x_1^0 и вычисляем значение функции в точке x ̅^=(x_1^0+∆x_1^0,x_2^0,…,x_n^0).
Если значение функции в данной точке больше, чем значение функции f(x ̅^0), то делаем пробный шаг по этой же переменной, но в противоположном направлении. Если значение функции и в точке x ̅^=(x_1^0-∆x_1^0,x_2^0,…,x_n^0) больше чем f(x ̅^0), то оставляем точку x ̅^0 без изменений. Иначе заменяем точку x ̅^0 на x ̅^ или наx ̅^ в зависимости от того, где значение функции меньше исходного. Из вновь полученной точки делаем пробные шаги по оставшимся координатам, используя тот же самый алгоритм.
Если в процессе исследующего поиска не удается сделать ни одного удачного пробного шага, то ∆x ̅ необходимо скорректировать (уменьшить). После чего вновь переходим к исследующему поиску.
Если в процессе исследующего поиска сделан хотя бы один удачный пробный шаг, то переходим к поиску по образцу.
Поиск по образцу. После исследующего поиска мы получаем точку x ̅^01. Направление x ̅^01-x ̅^0 определяет направление, в котором функция уменьшается. Поэтому проводим минимизацию функции в указанном на- правлении, решая задачу:
min┬λ〖f(x ̅^0+λ(x ̅^01-x ̅^0 ))〗.
В поиске по "образцу" величина шага по каждой переменной пропорциональна величине шага на этапе исследующего поиска. Если удастся сделать удачный шаг в поиске по "образцу", то в результате находим новое приближение x ̅^1=x ̅^0+λ(x ̅^01-x ̅^0 ), где
λ^0=arg min┬λ〖f(x ̅^0+λ(x ̅^01-x ̅^0 ))〗.
Из точки начинаем новый исследующий поиск и т.д.
Пример траектории спуска в алгоритме Хука-Дживса показан на рисунке 2.
Рисунок 2 – Пример траектории спуска в алгоритме Хука-Дживса
Метод Пауэлла
Метод Пауэлла относится к последовательным стратегиям. Задается начальная точка и с помощью пробного шага находится три точки так, чтобы они были как можно ближе к искомой точке минимума. В полученных точках вычисляются значения функции. Затем строится интерполяционный полином второй степени, проходящий через эти точки. В качестве приближения точки минимума берется точка минимума полинома. Поиск заканчивается, когда полученная точка отличается от лучшей из трех опорных точек не более чем на заданную величину [3].
Схему алгоритма можно описать следующим образом.
Пусть x_1 – начальная точка, х – выбранная величина шага по оси x.
Последовательность действий алгоритма Пауэлла следующая:
Вычислить x2 = x1 + ∆x.
Вычислить f(x1) и f(x2).
Если f(x1)> f(x2), положить x3 = x1 +2 ∆x (смотрите рисунок 3), иначе положить x3 = x1 - ∆x (смотрите рисунок 4).
Рисунок 3 – Отрезок xϵ[x1;x3]
Рисунок 4 – Отрезок xϵ[x3;x2]
Вычислить f(x3) и найти:
f_min=min{f_1,f_2,f_3 },
x_min – точка, которой соответствует f_min.
По трем точкам x_1, x_2 и x_3 вычислить x ̅, используя формулу квадратичной аппроксимации:
x ̅=1/2*(f(x_2^2-x_3^2 )+f(x_3^2-x_1^2 )+f(x_1^2-x_2^2 ))/(f(x_2-x_3 )+f(x_3-x_1 )+f(x_1-x_2 ) )
Также найти значение функции f(x ̅).
Если знаменатель в формуле для x ̅ на некоторой итерации обращается в нуль, то в этом случае результатом интерполяции является прямая, тогда следует положить x_1=x_min, и перейти к шагу 1.
Проверка на окончание поиска:
Являются ли разности f_min-f(x ̅) и x_min-x ̅ достаточно малыми?
Если оба условия выполняются, то необходимо закончить поиск, иначе переход к п.7.
Выбрать из x_min и x ̅ точку, в которой значение функции меньше, чем в другой, а также 2 точки по обе стороны от нее. Обозначить эти точки в естественном порядке и перейти к п.4.
Описание программы
Программа написана на языке высокого уровня Delphi (среда разработки: Delphi and C++ Builder 2009). Ее условно можно разделить на 3 части. Первая часть представляет собой графическое оформление программы и включает в себя прорисовку координатной сетки, а также иллюстрацию работы двух других частей программы как прорисовку пути движения исходной точки к точке минимума функции. Во второй части программы реализован метод Хука-Дживса с использованием метода Пауэлла. В третьей части реализован градиентный метод наискорейшего спуска с использованием метода Фибоначчи.
В данном разделе описана вторая часть программы.
Программа для вычисления точки минимума содержит в себе 2 основных метода: основную процедуру, в которой реализован метод Хука-Дживса, и вспомогательную функцию, в которой реализован метод Пауэлла.
Входными параметрами для реализации метода Хука-Дживса являются:
x0HJ и y0HJ – координаты начальной точки, относительно которой определяется начальное направление движения к точке минимума (вводится в интерфейсе программы).
delta – начальный шаг метода (вводится в интерфейсе программы).
EpsHJ – значение, до которого необходимо уменьшать delta для нахождения точки минимума.
В начале программы производится исследующий поиск.
Для этого ищется направление движения к точке минимума путем сравнения значений в восьми точках, радиально расположенных вокруг начальной точки на расстоянии delta, а также значения в самой начальной точке, которая находится в центре остальных восьми точек. В случае, если значение функции в начальной точке оказалось меньше, чем во всех остальных, значение переменной delta уменьшается в 2 раза и процесс поиска минимального значения функции в девяти точках повторяется до тех пор, пока delta не станет меньше EpsHJ.
Если же минимальное значение функции не в центральной точке, то обозначается направление движения поиска минимального значения функции путем присваивания переменным d1 и d2 значений 0,-1 или +1 в зависимости от направления движения по осям Х и У соответственно.
Далее производится поиск минимального значения функции на обозначенной прямой с помощью метода Пауэлла, то есть начинается поиск по образцу.
Входными данными для метода Пауэлла являются значение координат текущей точки и направление прямой для исследования. Выходными данными являются значения координат точки с минимальным значением на обозначенной прямой.
Вычисление координат выходной точки метода Пауэлла производится путем использования аппроксимирующей функции. Причем, вычисления производятся до тех пор, пока погрешность вычисления значения функции в выходной точке не будет меньше EpsP1, а погрешность вычисления координат выходной точки – меньше EpsP2.
Значения следующих переменных записываются в файл «Хука-Дживса.txt» на каждой итерации метода Хука-Дживса:
k – номер итерации;
x – первая координата очередной точки метода;
y – вторая координата очередной точки метода;
f – значение функции в точке (x, y);
d – значение переменной delta.
Результат выполнения программы
Проверим работу программы при начальных условиях x0=2, y0=-3 при разных значениях delta (смотрите рисунки 5-7).
Рисунок 5 – Результат выполнения программы при начальных условиях х0=2, у0=-3 и delta0=1
Содержание выходного файла «Хука-Дживса.txt»:
k x y f d
0 2.0000000000 -3.0000000000 33.0000000000 1.000
1 -0.7000000000 -0.3000000000 -3.4500000000 1.000
2 -0.7000000000 -0.3000000000 -3.4500000000 0.500
3 -1.0600000000 -0.6600000000 -4.0980000000 0.500
4 -1.0600000000 -0.6600000000 -4.0980000000 0.250
5 -1.0600000000 -0.6600000000 -4.0980000000 0.125
6 -1.0600000000 -0.7500000000 -4.1142000000 0.125
7 -1.0600000000 -0.7500000000 -4.1142000000 0.063
8 -1.0000000000 -0.7500000000 -4.1250000000 0.063
9 -1.0000000000 -0.7500000000 -4.1250000000 0.031
10 -1.0000000000 -0.7500000000 -4.1250000000 0.016
11 -1.0000000000 -0.7500000000 -4.1250000000 0.008
12 -1.0000000000 -0.7500000000 -4.1250000000 0.004
13 -1.0000000000 -0.7500000000 -4.1250000000 0.002
14 -1.0000000000 -0.7500000000 -4.1250000000 0.001
Рисунок 6 – Результат выполнения программы при начальных условиях х0=2, у0=-3 и delta0=0.01
Содержание выходного файла «Хука-Дживса.txt»:
k x y f d
0 2.0000000000 -3.0000000000 33.0000000000 0.010
1 -0.7000000000 -0.3000000000 -3.4500000000 0.010
2 -1.0600000000 -0.6600000000 -4.0980000000 0.010
3 -0.9880000000 -0.7320000000 -4.1239200000 0.010
4 -1.0024000000 -0.7464000000 -4.1249568000 0.010
5 -1.0024000000 -0.7464000000 -4.1249568000 0.005
6 -1.0024000000 -0.7500000000 -4.1249827200 0.005
7 -1.0024000000 -0.7500000000 -4.1249827200 0.003
8 -1.0000000000 -0.7500000000 -4.1250000000 0.003
9 -1.0000000000 -0.7500000000 -4.1250000000 0.001
10 -1.0000000000 -0.7500000000 -4.1250000000 0.001
Рисунок 7 – Результат выполнения программы при начальных условиях х0=2, у0=-3 и delta0=10
Содержание выходного файла «Хука-Дживса.txt»:
k x y f d
0 2.0000000000 -3.0000000000 33.0000000000 10.000
1 2.0000000000 -3.0000000000 33.0000000000 5.000
2 -1.0000000000 -3.0000000000 6.0000000000 5.000
3 -1.0000000000 -3.0000000000 6.0000000000 2.500
4 -1.0000000000 -0.7500000000 -4.1250000000 2.500
5 -1.0000000000 -0.7500000000 -4.1250000000 1.250
6 -1.0000000000 -0.7500000000 -4.1250000000 0.625
7 -1.0000000000 -0.7500000000 -4.1250000000 0.313
8 -1.0000000000 -0.7500000000 -4.1250000000 0.156
9 -1.0000000000 -0.7500000000 -4.1250000000 0.078
10 -1.0000000000 -0.7500000000 -4.1250000000 0.039
11 -1.0000000000 -0.7500000000 -4.1250000000 0.020
12 -1.0000000000 -0.7500000000 -4.1250000000 0.010
13 -1.0000000000 -0.7500000000 -4.1250000000 0.005
14 -1.0000000000 -0.7500000000 -4.1250000000 0.002
15 -1.0000000000 -0.7500000000 -4.1250000000 0.001
16 -1.0000000000 -0.7500000000 -4.1250000000 0.001
Как видно из рисунков 5-7 на быстроту сходимости сильно влияет величина delta. В данном случае самая быстрая сходимость получается при delta = 0.001 (10 итераций), самая медленная – при delta = 10 (16 итераций).
Метод наискорейшего спуска с использованием метода Фибоначчи
Метод наискорейшего спуска
Градиент функции в любой точке показывает направление наибольшего локального увеличения f(x ̅). Поэтому при поиске минимума f(x ̅), следует двигаться в направлении противоположном направлению градиента ∇f(x ̅) в данной точке, то есть в направлении наискорейшего спуска [1].
Итерационная формула процесса наискорейшего спуска имеет вид:
x ̅^(k+1)=x ̅^k-λ^k*∇f(x ̅^k)
Очевидно, что в зависимости от выбора параметра λ траектории спуска будут существенно различаться (смотрите рисунок ?????????). При большом значении λ траектория спуска будет представлять собой колебательный процесс, а при слишком больших λ процесс может расходиться. При выборе малых λ траектория спуска будет плавной, но и процесс будет сходиться очень медленно. Обычно выбирают из условия:
λ^k=argmin┬λ〖f(x ̅^k-λ^k*∇f(x ̅^k ) ),〗
решая одномерную задачу минимизации с использованием некоторого метода. В этом случае получаем алгоритм наискорейшего спуска.
Рисунок ???? – Траектории спуска в зависимости от величины шага
Если λ определяется в результате одномерной минимизации, то градиент в точке очередного приближения будет ортогонален направлению предыдущего спуска:
∇f(x ̅^(k+1))⊥∇f(x ̅^k)
Процесс наискорейшего спуска обычно быстро сходится вдали от точки экстремума и медленно в районе экстремума. Поэтому метод наискорейшего спуска нередко используют в комбинации с другими алгоритмами.
Метод Фибоначчи
Последовательность чисел Фибоначчи подчиняется соотношению [1]:
F_(n+2)=F_(n+1)+F_n,
где n=1,2,3,… и F1=F2=1.
Она имеет вид: 1, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233, 377, 610, 987, 1597 … .
С помощью индукции можно показать, что n-е число Фибоначчи вычисляется по формуле Бинэ:
F_n=(((1+√5)/2)^n-((1-√5)/2)^n)/√5,
где n=1,2,3,…
Из данной формулы видно, что при больших n значениях выполняется соотношение:
F_n≈((1+√5)/2)^n/√5,
так, что числа Фибоначчи с увеличением n растут очень быстро.
Сам алгоритм метода Фибоначчи очень похож на алгоритм метода золотого сечения. На начальном интервале вычисляются точки по следующим формулам:
x_1=a_0+F_n/F_(n+2) (b_0-a_0 ),
x_2=a_0+F_(n+1)/F_(n+2) (b_0-a_0 )=a_0+b_0-x_1.
Интервал неопределенности сокращается точно так же, как в методе золотого сечения. И на новой итерации вычисляется только одна новая точка и значение функции в ней.
На k-й итерации получаем точку минимальным значением, которая совпадает с одной из точек, вычисляемых по формулам:
x_1^=a_k+F_(n-k+1)/F_(n-k+3) (b_k-a_k ),
x_2^=a_k+F_(n-k+2)/F_(n-k+3) (b_k-a_k ).
и расположены на отрезке [ak, bk] симметрично относительно его середины.
Нетрудно заметить, что при k=n точки
x_1^=a_n+F_1/F_(n+2) (b_0-a_0 ) и x_2^=a_n+F_2/F_(n+2) (b_0-a_0 )
совпадают и делят отрезок [an, bn] пополам.
Следовательно,
(b_n-a_n)/2=(b_n-a_n)/F_(n+2) <.
Отсюда можно выбрать n из условия:
(b_n-a_n)/<F_(n+2).
Таким образом, это условие позволяет до начала работы алгоритма определить число итераций, необходимое для определения минимума с точностью при начальной величине интервала [a0, b0].
С ростом n из-за того, что Fn/Fn+2 – бесконечная десятичная дробь, возможно “искажение” метода и потеря интервала с точкой минимума (вследствие погрешностей вычислений).
Определение отрезка унимодальности
Функция называется унимодальной на некотором отрезке [a, b], если она имеет только одну точку минимума х* на отрезке [a, b], причем функция строго убывает на полуинтервале [a, х*) и строго возрастает на полуинтервале (х*, b]. Функция на отрезке [a, b] не обязательно должна быть непрерывной и дифференцируемой [2].
Для определения отрезка унимодальности будем использовать метод Свенна.
Перед вычислениями необходимо задать начальную точку х0, а также начальный шаг t>0.
Алгоритм метода Свенна следующий:
Вычислить значения заданной функции f(x) в точках х0-t, х0 и х0+t.
Если f(х0-t)>f(х0)<f(х0+t), то отрезок унимодальности найден. Если f(х0-t)<f(х0)>f(х0+t), то отрезок не унимодальный. Если не выполняется ни одно из условий из п.2, перейти к п.3.
Если f(х0-t)>f(х0)>f(х0+t), то функция убывающая, определить переменные: a=x0, ∆=t, x1=x0+∆, k=1. Если f(х0-t)<f(х0)<f(х0+t), то функция возрастающая, определить переменные: b=x0, ∆=-t, x1=x0+∆,k=1.
Вычислить следующую точку по формуле: xk+1=xk+2k ∆.
Если ∆>0 и f(xk)>f(xk+1), то определить переменные: a=xk, k=k+1, перейти к п.4. Если ∆<0 и f(xk)>f(xk+1), то определить переменные: b=xk, k=k+1, перейти к п.4. Если не выполняется ни одно из условий из п.5, то перейти к п.6.
Если ∆>0, то b= xk+1, иначе a= xk+1. Отрезок унимодальности найден.
Описание программы
Метод наискорейшего спуска с использованием метода Фибоначчи реализован в третьей части программы (смотрите п.4.3).
Эта часть программы состоит из трех основных методов: процедуры, реализующей метод наискорейшего спуска, функции, реализующей метод Фибоначчи, и функции, реализующей метод Свенна.
Входными данными для метода наискорейшего спуска являются:
x0F и y0F – координаты начальной точки (вводится в интерфейсе программы);
EpsF1 – погрешность вычисления значения функции в точке минимума;
EpsF2 – погрешность вычисления координат точки минимума.
На каждой итерации производится поиск направления движения к точке минимума. Это направление ищется с помощью функции, которая считает значение производной функции в указанной точке. Выражение производной функции задается явно в программе, поэтому при изменении выражения функции в программе также необходимо изменять выражение производной функции. После нахождения направления поиска точки минимума производится поиск этой точки с помощью метода Фибоначчи.
Итерации продолжаются до тех пор, пока расстояние между текущей и предыдущей точками минимума не будет меньше значения EpsF1, а значения функций в этих точках будут отличаться меньше, чем на EpsF2.
Функция, реализующая метод Фибоначчи в качестве входных параметров требует координаты предыдущей точки, а возвращает координаты текущей точки минимума.
В начале работы процедуры производится поиск отрезка унимодальности с помощью функции, реализующей метод Свенна. После этого производится поиск точки минимума методом Фибоначчи.
Функция, реализующая метод Свенна, в качестве входных параметров требует значения координат начальной точки, а также начальный шаг вычислений. На выходе же функция возвращает отрезок унимодальности в виде структуры, состоящей из массива двух точек.
На каждой итерации метода наискорейшего спуска производится запись следующих значений переменных в файл «Фибоначчи.txt»:
k – номер итерации;
x – первая координата очередной точки метода;
y – вторая координата очередной точки метода;
f – значение функции в точке (x, y).
Описание работы программы
Проверим работу программы при начальных условиях x0=2, y0=-3 (смотрите рисунок ??????).
Рисунок ??? – Результат выполнения программы при начальных условиях х0=2 и у0=-3
Содержание выходного файла «Фибоначчи.txt»:
k x y f
0 2.0000000000 -3.0000000000 33.0000000000
1 -1.2348152786 -1.3825923607 -3.1592391652
2 -0.9216302170 -0.8201119899 -4.0967431491
3 -1.0067230183 -0.7693609189 -4.1241147127
4 -0.9976551527 -0.7519518522 -4.1249758856
5 -1.0001835353 -0.7505487675 -4.1249992967
6 -0.9999346076 -0.7500525745 -4.1249999816
7 -1.0000046270 -0.7500150448 -4.1249999995
Выводы
В ходе выполнения данной курсовой работы были изучены основные методы оптимизации функции нескольких переменных, а также написана программа для вычисления точки минимума заданной функции.
Кроме того, заданная функция была решена аналитически и с использованием разработанной программы. Результаты вычислений на программе совпали с результатами, полученными решением задачи аналитическим методом. Этот факт свидетельствует о том, что алгоритмы программы были реализованы верно и программа работает правильно.
Сравнительный анализ двух алгоритмов, реализованных в разработанной программе показывает, что при вычислении градиентным методом наискорейшего спуска нахождение оптимальной точки производится за меньшее количество итераций, чем при вычислении прямым методом Хука-Дживса, в среднем в 2 раза быстрее.
Также было установлено, что на количество итераций в обоих методах влияет положение начальной точки, а также допустимые погрешности измерений.
Приложение. Программа на ЯВУ
unit Unit1;
interface
uses
Windows, Messages, SysUtils, Variants, Classes, Graphics, Controls, Forms,
Dialogs, StdCtrls, ExtCtrls;
const
k0 = -3;
indCL = 1;
countCL = 5;
EpsCL = 0.05;
EpsHJ = 0.001;
EpsP1 = 0.001;
EpsP2 = 0.001;
EpsF = 0.001;
EpsF1 = 0.0001;
EpsF2 = 0.0001;
FibNumsCount = 31;
type
TForm1 = class(TForm)
Panel1: TPanel;
Image1: TImage;
Button1: TButton;
Button2: TButton;
Button3: TButton;
Label1: TLabel;
Label2: TLabel;
Button4: TButton;
Label3: TLabel;
LabeledEdit1: TLabeledEdit;
LabeledEdit2: TLabeledEdit;
procedure Button1Click(Sender: TObject);
procedure Button2Click(Sender: TObject);
procedure Button3Click(Sender: TObject);
procedure Button4Click(Sender: TObject);
private
{ Private declarations }
public
{ Public declarations }
end;
type Point = record //структура, состоящая
x: Real; //из координат точки
y: Real;
f: Real; //и значения функции в этой точке
end;
type Segment = Array[1..2] of Real; //вершины отрезка
var
Form1: TForm1;
x0: Integer; //точка начала координат
y0: Integer;
iX: Integer; //вспомогательные переменные
iY: Integer;
kiY: Real;
kiX: Real;
NinePoints: Array[1..3,1..3] of Point; //для метода Хука-Дживса
FibNums: Array[1..FibNumsCount] of LongInt; //для метода Фибоначчи
implementation
{$R *.dfm}
function func(x: Real; y:Real): Real;
//высчитывание значения функции в координатах х и у
begin
func := 3*x*x+2*y*y+6*x+3*y;
end;
function derX(x: Real; y: Real): Real;
//высчитывание значения производной функции в координатах х и у
//по координате х
begin
derX:=6*x+6;
end;
function derY(x: Real; y: Real): Real;
//высчитывание значения производной функции в координатах х и у
//по координате у
begin
derY:=4*y+3;
end;
function funcP(p: Point): Real;
//высчитывание значения функции в координатах р.х и р.у
begin
funcP := func(p.x,p.y);
end;
//4 функции для рисования
function J_To_Y(j: Integer): Real;
begin
J_To_Y:=(y0-j)/kiY;
end;
function I_To_X(i: Integer): Real;
begin
I_To_X:=(i-x0)/kiX;
end;
function X_To_I(x: Real): Integer;
begin
X_To_I:=round(x*kiX+x0);
end;
function Y_To_J(y: Real): Integer;
begin
Y_To_J:=round(y0-y*kiY);
end;
procedure TForm1.Button1Click(Sender: TObject);
//прорисовка сетки
var
i:Integer;
begin
with Image1 do
begin
iY := trunc(abs(StrToFloat(LabeledEdit1.Text))+1);
iX := trunc(abs(StrToFloat(LabeledEdit2.Text))+1);
if iX<3 then
iX:=3;
if iX>iY then
iY := iX
else
iX := iY;
x0 := Width div 2;
y0 := Height div 2;
kiY := y0/iY;
kiX := x0/iX;
with Canvas do
begin
Brush.Color:=clWhite;
Pen.Color:=clWhite;
Rectangle(0,0,Image1.Width,Image1.Height);
//рисование осей
Pen.Color:=clBlack;
Pen.Width:=2;
MoveTo(Image1.Left,y0);
LineTo(Image1.Width,y0);
MoveTo(x0,Image1.Top);
LineTo(x0,Image1.Height);
//штрихи и подписи на осях
Pen.Color:=$00ff00;
Pen.Width:=1;
for i := 1 to iY-1 do
begin
MoveTo(0,y0-round(kiY*i));
LineTo(Image1.Width,y0-round(kiY*i));
TextOut(x0-15,y0-round(kiY*i)-15,IntToStr(i));
MoveTo(0,y0+round(kiY*i));
LineTo(Image1.Width,y0+round(kiY*i));
TextOut(x0-15,y0+round(kiY*i)-15,IntToStr(-i));
end;
for i := 1 to iX-1 do
begin
MoveTo(x0-round(kiX*i),0);
LineTo(x0-round(kiX*i),Image1.Height);
TextOut(x0-round(kiX*i)-15,y0-15,IntToStr(-i));
MoveTo(x0+round(kiX*i),0);
LineTo(x0+round(kiX*i),Image1.Height);
TextOut(x0+round(kiX*i)-15,y0-15,IntToStr(i));
end;
end;
end;
end;
procedure TForm1.Button2Click(Sender: TObject);
//линии уровня
var
i: Integer;
j: Integer;
k: Integer;
x: Real;
y: Real;
begin
Button1Click(Self);
with Image1 do
begin
Canvas.Pen.Color:=clGray;
Canvas.Pen.Width:=2;
k:=k0;
while k-k0<countCL*indCL do
begin
for i := 0 to Width do
begin
x:=I_To_X(i);
for j := 0 to Height do
begin
y:=J_To_Y(j);
if (abs(func(x,y)-k)<EpsCL) then
begin
Canvas.MoveTo(i,j);
Canvas.LineTo(i,j);
end;
end;
end;
k:=k+indCL;
end;
end;
end;
function min(p: Point; delta: Real):Point;
//нахождение точки с наименьшим значением функции среди девяти
var
i: Integer;
j: Integer;
minim: Point;
begin
for i := 1 to 3 do
for j := 1 to 3 do
with NinePoints[i,j] do
begin
x:=p.x+(i-2)*delta;
y:=p.y-(j-2)*delta;
f:=func(x,y);
end;
minim:=NinePoints[1,1];
for i := 1 to 3 do
for j := 1 to 3 do
if NinePoints[i,j].f<minim.f then
minim:=NinePoints[i,j];
min:=minim;
end;
function Pauella(pj: Point; d1: Integer; d2: Integer): Point;
//Метод Пауэлла
var
i: Integer;
imin: Integer;
lmin: Real;
fmin: Real;
delt: Real;
l0: Real;
lW: Real;
fW: Real;
flg: Boolean;
flg1: Boolean;
l: Array[1..3] of Real;
f: Array[1..3] of Real;
function funcL(l:Real):Real;
//нахождение значения функции на обозначенной прямой
begin
funcL:=func(pj.x+l*d1, pj.y+l*d2);
end;
function minF:Integer;
//нахождение точки с минимальным значением функции среди трех
var
minim: Real;
ind: Integer;
i:Integer;
begin
minim:=f[1];
ind:=1;
for i := 2 to 3 do
if f[i]<minim then
begin
minim:=f[i];
ind:=i;
end;
minF:=ind;
end;
function funcLW:Real;
//нахождение точки с минимальным значением на обозначенной прямой
begin
if f[1]*(l[2]-l[3])+f[2]*(l[3]-l[1])+f[3]*(l[1]-l[2])=0 then
begin
flg:=True;
l0:=l[minF];
end
else
funcLW:=(f[1]*(l[2]*l[2]-l[3]*l[3])+f[2]*(l[3]*l[3]-l[1]*l[1])+
f[3]*(l[1]*l[1]-l[2]*l[2]))/
(2*(f[1]*(l[2]-l[3])+f[2]*(l[3]-l[1])+f[3]*(l[1]-l[2])));
end;
begin
delt:=1;
l0:=0;
flg:=True;
repeat
flg1:=False;
if flg then
begin
l[1]:=l0;
l[2]:=l[1]+delt;
for i := 1 to 2 do
f[i]:=funcL(l[i]);
if f[1]>f[2] then
begin
l[3]:=l[1]+2*delt;
f[3]:=funcL(l[3]);
end
else
begin
for i := 3 downto 2 do
begin
f[i]:=f[i-1];
l[i]:=l[i-1];
end;
l[1]:=l[2]-delt;
f[1]:=funcL(l[1]);
end;
end;
imin:=minF;
fmin:=f[imin];
lmin:=l[imin];
lW:=funcLW;
if flg1 then
continue;
fW:=funcL(lW);
if (lW>=l[1])and(lW<=l[3]) then
begin
flg:=False;
if lW<l[2] then
if fW<fmin then
begin
l[3]:=l[2];
l[2]:=lW;
end
else
l[1]:=lW
else
if fW<fmin then
begin
l[1]:=l[2];
l[2]:=lW;
end
else
l[3]:=lW;
for i := 1 to 3 do
f[i]:=funcL(l[i]);
end
else
begin
l0:=lW;
flg:=True;
end;
until ((abs(fmin-fW)<EpsP1) and (abs(lmin-lW)<EpsP2));
Pauella.x:=pj.x+lW*d1;
Pauella.y:=pj.y+lW*d2;
end;
procedure TForm1.Button3Click(Sender: TObject);
//Хука-Дживса
var
d1: Integer;
d2: Integer;
delta: Real;
x0HJ: Real; //координаты
y0HJ: Real; //начальной точки
p: Point; //текущая точка
pj: Point; //предыдущая точка
kit: Integer;
flag: Boolean;
f: TextFile;
begin
AssignFile(f,Хука-Дживса.txt);
Rewrite(f);
Writeln(f, k:4,x:15,y:15,f:15,d:7);
x0HJ := StrToFloat(LabeledEdit1.Text);
y0HJ := StrToFloat(LabeledEdit2.Text);;
Image1.Canvas.Pen.Color:=clRed;
Image1.Canvas.Pen.Width:=2;
kit:=0;
delta:=1;
flag:=False;
pj.x:=x0HJ;
pj.y:=y0HJ;
Image1.Canvas.MoveTo(X_To_I(pj.x),Y_To_J(pj.y));
pj.f:=funcP(pj);
Writeln(f,kit:4,pj.x:15:10,pj.y:15:10,pj.f:15:10,delta:7:3);
while not flag do
begin
inc(kit);
p:=min(pj,delta);
if ((p.x=pj.x) and (p.y=pj.y)) then
begin
delta:=delta/2;
Writeln(f,kit:4,pj.x:15:10,pj.y:15:10,funcP(pj):15:10,delta:7:3);
if delta<EpsHJ then
flag:=True;
continue;
end;
d1:=round((p.x-pj.x)/abs(p.x-pj.x+0.001));
d2:=round((p.y-pj.y)/abs(p.y-pj.y+0.001));
pj:=Pauella(pj,d1,d2);
Writeln(f,kit:4,pj.x:15:10,pj.y:15:10,funcP(pj):15:10,delta:7:3);
Image1.Canvas.LineTo(X_To_I(pj.x),Y_To_J(pj.y));
end;
Label1.Caption:=Метод Хука-Дживса: x=+
FloatToStr(pj.x)+
; y=+FloatToStr(pj.y)+
; f=+FloatToStr(funcP(pj))+
; количество итераций=+IntToStr(kit);
CloseFile(f);
end;
function funcF(x:Real; y:Real; l:Real):Real;
begin
funcF:=func(x-l*derX(x,y),y-l*derY(x,y));
end;
function Svenna(p0:Point; step:Real):Segment;
//Свенна
function funcP0(l:Real): Real;
//высчитывание значения функции в точке р0
begin
funcP0:=funcF(p0.x,p0.y,l);
end;
var
i: Integer;
delt: Real;
l: Array[1..3] of Real;
f: Array[1..3] of Real;
begin
delt:=step;
l[1]:=-abs(delt);
l[2]:=0;
l[3]:=abs(delt);
for i := 1 to 3 do
f[i]:=funcP0(l[i]);
if (f[1]<f[2]) or (f[2]>f[3]) then
if (f[1]<f[2]) and (f[2]>f[3]) then
begin
Form1.Label1.Caption:=Отрезка унимодальности не существует;
exit;
end
else
repeat
if f[2]>f[3] then
begin
for i := 1 to 2 do
begin
l[i]:=l[i+1];
f[i]:=f[i+1];
end;
delt:=abs(delt*2);
l[3]:=l[3]+delt;
f[3]:=funcP0(l[3]);
end
else
begin
for i := 3 downto 2 do
begin
l[i]:=l[i-1];
f[i]:=f[i-1];
end;
delt:=-abs(delt*2);
l[1]:=l[1]-delt;
f[1]:=funcP0(l[1]);
end;
until (f[1]>f[2]) and (f[2]<f[3]);
Svenna[1]:=l[1];
Svenna[2]:=l[3];
end;
function Fibonacci(p:Point):Point;
//Метод Фибонначи
function funcP0(l:Real): Real;
//высчитывание значения функции в точке р0
begin
funcP0:=funcF(p.x,p.y,l);
end;
var
i: Integer;
indF: Integer;
lSeg: Segment;
l: Array[1..2] of Real;
f: Array[1..2] of Real;
begin
lSeg:=Svenna(p,1);
indF:=1;
while FibNums[indF]<(abs(lSeg[2]-lSeg[1]/EpsF)) do
inc(indF);
l[1]:=lSeg[1]+(lSeg[2]-lSeg[1])*FibNums[indF-2]/FibNums[indF];
l[2]:=lSeg[1]+(lSeg[2]-lSeg[1])*FibNums[indF-1]/FibNums[indF];
repeat
for i := 1 to 2 do
f[i]:=funcP0(l[i]);
if f[1]>f[2] then
begin
lSeg[1]:=l[1];
l[1]:=l[2];
l[2]:=lSeg[1]+lSeg[2]-l[1];
end
else
begin
lSeg[2]:=l[2];
l[2]:=l[1];
l[1]:=lSeg[1]+lSeg[2]-l[2];
end;
until abs(lSeg[2]-lSeg[1])<EpsF;
lSeg[1]:=lSeg[1]+(lSeg[2]-lSeg[1])/2;
Fibonacci.x:=p.x-lSeg[1]*derX(p.x,p.y);
Fibonacci.y:=p.y-lSeg[1]*derY(p.x,p.y);
end;
procedure TForm1.Button4Click(Sender: TObject);
//Метод наискорейшего спуска
var
i: Integer;
kit: Integer;
x0F: Real;
y0F: Real;
f: TextFile;
p: Array[1..2] of Point;
begin
AssignFile(f,Фибоначчи.txt);
Rewrite(f);
Writeln(f, k:4,x:15,y:15,f:15);
x0F := StrToFloat(LabeledEdit1.Text);
y0F := StrToFloat(LabeledEdit2.Text);
Image1.Canvas.Pen.Color:=clBlue;
Image1.Canvas.Pen.Width:=3;
kit:=0;
p[2].x:=x0F;
p[2].y:=y0F;
p[2].f:=funcP(p[2]);
Writeln(f,kit:4,p[2].x:15:10,p[2].y:15:10,p[2].f:15:10);
Image1.Canvas.MoveTo(X_To_I(p[2].x),Y_To_J(p[2].y));
FibNums[1]:=1;
FibNums[2]:=1;
for i := 3 to FibNumsCount do
FibNums[i]:=FibNums[i-2]+FibNums[i-1];
repeat
inc(kit);
p[1]:=p[2];
p[2]:=Fibonacci(p[1]);
Writeln(f,kit:4,p[2].x:15:10,p[2].y:15:10,funcP(p[2]):15:10);
Image1.Canvas.LineTo(X_To_I(p[2].x),Y_To_J(p[2].y));
until ((Sqrt(Sqr(p[2].x-p[1].x)+Sqr(p[2].y-p[1].y))<EpsF1) and
(abs(p[2].f-p[1].f)<EpsF2));
Label2.Caption:=Метод Фибоначчи: x=+
FloatToStr(p[2].x)+ ; y=+FloatToStr(p[2].y)+
; f=+FloatToStr(funcP(p[2]))+
; количество итераций=+IntToStr(kit);
CloseFile(f);
end;
end.
Список используемых источников
Б.Ю. Лемешко «Методы оптимизации: Конспект лекций». – Новосибирск: Изд-во НГТУ, 2009. – 126 с.
И.Н. Мастяева, О.Н. Семенихина «Методы оптимизации» ¬– МЭСИ, 2003. – 135 с.
Методы Оптимизации Систем Автоматизированного Проектирования. URL: http://optimizaciya-sapr.narod.ru/bez_odnomer/paul.html (дата просмотра: 31 марта 2014 г.)