Оптимизация методом покоординатного спуска (Гаусса-Зейделя)

@Сергей.MS · Регистрация: 16.03.2012

Студворк — интернет-сервис помощи студентам

Есть рабочий вариант:

Matlab M

clear all; clc;
% Значения коэффициентов
c1 = -2; 
c2 = -1; 
c12 = 1; 
g = 0.2; % постоянная шага
d = 0.01; % дельта
% Начальная точка
x1 = 16; 
x2 = 18; 
k = 1;  % Счетчик шагов
kmax = 100; % Предельное число шагов, 
% задается для предотвращения зацикливания
% Массивы для хранения промежуточных координат
x1trace = [x1]; 
x2trace = [x2]; 
i = 2; 
while k < kmax 
% Спуск по первой координате
gr1 = 2*c1*x1 + c12*x2; 
x1 = x1 + g*gr1; 
% Сохранение координат
x1trace(i) = x1; 
x2trace(i) = x2; 
i = i + 1; 
% Спуск по второй координате
gr2 = 2*c2*x2 + c12*x1; 
x2 = x2 + g*gr2; 
% Сохранение координат
x1trace(i) = x1; 
x2trace(i) = x2; 
i = i + 1; 
% Проверка условия останова
if sqrt(gr1^2 + gr2^2) <= d; 
break; % Выход из цикла в случае выполнения условия
end
k = k + 1; 
end
% Построение графика
x = -20:0.1:20; 
y = -20:0.1:20; 
[X, Y] = meshgrid(x, y); 
Z = c1*X.^2 + c2*Y.^2 + c12*X.*Y; 
[C, h] = contour(X, Y, Z); 
clabel(C, h); 
% Отображение меток на линиях уровня
hold on; 
plot(x1trace, x2trace, '-'); 
% Вывод начальной точки на график
text(x1trace(1) + 0.2, x2trace(1) + 0.5, 'M0'); 
% Вывод решения на график
text(x1 + 2, x2, ...
strvcat(['x1 = ' (num2str(x1))], ...
        ['x2 = ' (num2str(x2))], ...
        ['k = '  (num2str(k))]));

Необходимо по аналогии построить график для другой функции:
Начал изменять код под свою функцию, пока не получается

Matlab M

clear all; clc;
% Значения коэффициентов
c0  = -0.8; 
c1  =  1.3; 
c11 = -6.1;
c12 = -0.6;
c2  =  1.7;
c22 = -1.7;
g = 0.2; % постоянная шага
d = 0.01; % дельта
% Начальная точка
x1 = 16; 
x2 = 18; 
k = 1;  % Счетчик шагов
kmax = 100; % Предельное число шагов, 
% задается для предотвращения зацикливания
% Массивы для хранения промежуточных координат
x1trace = [x1]; 
x2trace = [x2]; 
i = 2; 
while k < kmax 
% Спуск по первой координате
gr1 = c0 + c1*x1 + c12*x1 + 2*c11*x1; 
x1 = x1 + g*gr1; 
% Сохранение координат
x1trace(i) = x1; 
x2trace(i) = x2; 
i = i + 1; 
% Спуск по второй координате
gr2 = c0 + c2*x2 + c12*x2 + 2*c22*x2; 
x2 = x2 + g*gr2; 
% Сохранение координат
x1trace(i) = x1; 
x2trace(i) = x2; 
i = i + 1; 
% Проверка условия останова
if sqrt(gr1^2 + gr2^2) <= d; 
break; % Выход из цикла в случае выполнения условия
end
k = k + 1; 
end
% Построение графика
x = -60:0.1:60; 
y = -60:0.1:60; 
[X, Y] = meshgrid(x, y); 
Z = c0 + c1*X + c2*Y + c12*X.*Y + c11*X.^2 + c22*Y.^2; 
[C, h] = contour(X, Y, Z); 
clabel(C, h); 
% Отображение меток на линиях уровня
hold on; 
plot(x1trace, x2trace, '-'); 
% Вывод начальной точки на график
text(x1trace(1) + 0.2, x2trace(1) + 0.5, 'M0'); 
% Вывод решения на график
text(x1 + 2, x2, ...
strvcat(['x1 = ' (num2str(x1))], ...
        ['x2 = ' (num2str(x2))], ...
        ['k = '  (num2str(k))]));

на рисунке 1 - рабочая функция
2 - необходимо разработать под эту функцию

@Сергей.MS · 13.01.2013, 23:42 **[ТС]**

по читал инфу:
итог - нужно в gr1 и gr2 прописать производные соответственно по x1 и x2 от функции f(x1,x2)
производные нашёл, прописал:
gr1 = c1 + c12*x2 + 2*c11*x1
gr2 = c2 + c12*x1 + 2*c22*x2
но всё равно пока ни как, график должен сходиться в точку, ну а у меня расходится (значения gr1 и gr2 должны уменьшаться)
если кто волокёт в этом деле помогите(((

@Сергей.MS · 18.01.2013, 14:09 **[ТС]**

Вот правильное решение поставленной задачи:

Matlab M

clear all; clc;
% Значения коэффициентов
c0  =  0.8; 
c1  =  1.3; 
c11 = -6.1;
c12 = -0.6;
c2  =  1.7;
c22 = -1.7;
g = 0.05; % постоянная шага
d = 0.01; % дельта
% Начальная точка
x1 = -9; 
x2 = 8; 
k = 1;  % Счетчик шагов
kmax = 100; % Предельное число шагов, 
% задается для предотвращения зацикливания
% Массивы для хранения промежуточных координат
x1trace = [x1]; 
x2trace = [x2]; 
i = 2; 
while k < kmax 
% Спуск по первой координате
gr1 = c1 + c12*x2 + 2*c11*x1; 
x1 = x1 + g*gr1
% Сохранение координат
x1trace(i) = x1; 
x2trace(i) = x2; 
i = i + 1; 
% Спуск по второй координате
gr2 = c2 + c12*x1 + 2*c22*x2; 
x2 = x2 + g*gr2
% Сохранение координат
x1trace(i) = x1; 
x2trace(i) = x2; 
i = i + 1; 
% Проверка условия останова
if sqrt(gr1^2 + gr2^2) <= d; 
break; % Выход из цикла в случае выполнения условия
end
k = k + 1; 
end
% Построение графика
x = -10:0.1:10; 
y = -10:0.1:10; 
[X, Y] = meshgrid(x, y); 
Z = c0 + c1*X + c2*Y + c12*X.*Y + c11*X.^2 + c22*Y.^2; 
[C, h] = contour(X, Y, Z); 
clabel(C, h); 
% Отображение меток на линиях уровня
hold on;
plot(x1trace, x2trace, '-'); 
% Вывод начальной точки на график
text(x1trace(1) + 0.2, x2trace(1) + 0.5, 'M0'); 
% Вывод решения на график
text(x1 + 2, x2, ...
strvcat(['x1 = ' (num2str(x1))], ...
        ['x2 = ' (num2str(x2))], ...
        ['k = '  (num2str(k))])); 
hold off;

В итоге получилось см ниже

Зосима · 18.01.2013, 14:35

Сергей.MS, супер! Особенно график!

А как этот метод называется?

@Сергей.MS · 18.01.2013, 14:37 **[ТС]**

Метод покоординатного спуска (Гаусса-Зейделя)

Зосима · 18.01.2013, 14:51

А в чем там загвоздка была? Как переписать ее для другой ф-ции?

@Сергей.MS · 18.01.2013, 14:56 **[ТС]**

Вообщем сходил на консультацию в универ, и вот что узнал
вся замануха вот в чем:
1. необходимо правильно найти производные от функции двух переменных gr1 b gr2, не сразу получилось)
2. подобрать самостоятельно постоянную шага, так чтобы x1'' было меньше x1' и x2'' меньше x2'.

P.S. зря переименовал тему, мне ещё 3 метода решить надо будет))

Зосима · 18.01.2013, 15:05

Спасибо!

Сообщение от Сергей.MS

мне ещё 3 метода решить надо будет))

Они будут каждый в своей теме!

@Сергей.MS · 18.01.2013, 21:59 **[ТС]**

по ходу возник вопрос, а как изменить код, чтобы график выводился на экран с одинаковым масштабом по осям, если одинаковые величины откладываются, см график. по x -+10 и по y -+10 а в итоге на экране прямоугольник?

Зосима · 18.01.2013, 22:43

после plot напиши
axes equal

@mavi93 · 22.05.2014, 00:06

Вот вы сделали с постоянным шагом, а по методу оптимизации Гаусса-Зейделя шаг этот зависит от предыдущих итераций, как это сделать

вот ссылочка теорию и пример

@olegar4 · 26.11.2014, 00:38

этот код работает только для определенной функции, так? А можно создать М-файл для общего решения?

Зосима · 26.11.2014, 12:40

Сообщение от olegar4

А можно создать М-файл для общего решения?

почему бы и нет

Matlab M

function Xm = GaussSeidelMin(fun, x0, d, g)
% минимизация покоординатным спуском (метод Гаусса-Зейделя)
% Xm - координтаы минимума
% fun - функция двух переменных, котой ищется миниммум
% x0 - координаты начальной точки
% d - точность
% g - постоянная шага
 
if nargin<=4
    g = 0.05; % значение по умолчанию
end
if nargin<=3
    d = 0.001; % значение по умолчанию
end
 
k = 1;  % Счетчик шагов
kmax = 100; % Предельное число шагов,
% задается для предотвращения зацикливания
% Массивы для хранения промежуточных координат
x1 = x0(1);
x2 = x0(2);
x1trace = [x1];
x2trace = [x2];
i = 2;
h = 1e-6; % приращение для рассчета градиентов
while k < kmax
    % Спуск по первой координате
    gr1 = ( fun(x1+h,x2)-fun(x1-h,x2) )/(2*h); % производная по х1
    x1 = x1 - g*gr1;
    % Сохранение координат
    x1trace(i) = x1;
    x2trace(i) = x2;
    i = i + 1;
    % Спуск по второй координате
    gr2 = ( fun(x1,x2+h)-fun(x1,x2-h) )/(2*h); % производная по х2
    x2 = x2 - g*gr2;
    % Сохранение координат
    x1trace(i) = x1;
    x2trace(i) = x2;
    i = i + 1;
    % Проверка условия останова
    if sqrt(gr1^2 + gr2^2) <= d;
        break; % Выход из цикла в случае выполнения условия
    end
    k = k + 1;
end
 
Xm = [x1, x2]; % координаты минимума
if k >= kmax
    disp('Минимум не найден!')
    disp('Измените значение шага')
end
 
% Построение графика
d1 = abs(x0(1)-x1);
d2 = abs(x0(2)-x2);
x = linspace(x1-d1,x1+d1,50);
y = linspace(x2-d2,x2+d2,50);
[X, Y] = meshgrid(x, y);
Z = fun(X,Y);
[C, h] = contour(X, Y, Z);
clabel(C, h);
% Отображение меток на линиях уровня
hold on;
plot(x1trace, x2trace, '-');
% Вывод начальной точки на график
text(x1trace(1) + 0.2, x2trace(1) + 0.5, 'M0');
% Вывод решения на график
text(x1 + 0.3*d1, x2, ...
    strvcat(['x1 = ' (num2str(x1))], ...
    ['x2 = ' (num2str(x2))], ...
    ['k = '  (num2str(k))]));
hold off
end

А вот пример использования:

Matlab M

>> f = @(x,y) -exp( -x.^2 - y.^2 - 2*x );
>>GaussSeidelMin(f, [1 1], 1e-3, 0.1)
 
ans =
 
   -0.9976    0.0009

Оптимизация методом покоординатного спуска (Гаусса-Зейделя)

@mr_Adrian · 12.12.2014, 02:57

Сообщение от Зосима

gr1 = ( fun(x1+h,x2)-fun(x1-h,x2) )/(2*h); % производная по х1

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

Зосима · 12.12.2014, 10:49

mr_Adrian, по определению производная это предел отношения приращения ф-ции к приращению аргумента. Но в данном случае я нахожу не предел, а отношение при очень малом приращении h

это простейшая формула численного расчета производной

@mr_Adrian · 15.12.2014, 01:03

спасибо!. признателен.

@droi · 20.12.2014, 16:43

Не подскажите стандартные функции оптимизации в matlab?

Зосима · 20.12.2014, 17:55

droi, это fminsearch и fminbnd

@droi · 20.12.2014, 20:24

Зосима, Спасибо!

@apreliya1991 · 13.03.2015, 14:25

Пример выдает ошибку:

Error using Untitled2 (line 20)
Not enough input arguments.

Matlab M

function Xm = GaussSeidelMin(fun, x0, d, g)
% минимизация покоординатным спуском (метод Гаусса-Зейделя)
% Xm - координтаы минимума
% fun - функция двух переменных, котой ищется миниммум
% x0 - координаты начальной точки
% d - точность
% g - постоянная шага
 
if nargin<=4
    g = 0.05; % значение по умолчанию
end
if nargin<=3
    d = 0.001; % значение по умолчанию
end
 
k = 1;  % Счетчик шагов
kmax = 100; % Предельное число шагов,
% задается для предотвращения зацикливания
% Массивы для хранения промежуточных координат
x1 = x0(1);
x2 = x0(2);
x1trace = [x1];
x2trace = [x2];
i = 2;
h = 1e-6; % приращение для рассчета градиентов
while k < kmax
    % Спуск по первой координате
    gr1 = ( fun(x1+h,x2)-fun(x1-h,x2) )/(2*h); % производная по х1
    x1 = x1 - g*gr1;
    % Сохранение координат
    x1trace(i) = x1;
    x2trace(i) = x2;
    i = i + 1;
    % Спуск по второй координате
    gr2 = ( fun(x1,x2+h)-fun(x1,x2-h) )/(2*h); % производная по х2
    x2 = x2 - g*gr2;
    % Сохранение координат
    x1trace(i) = x1;
    x2trace(i) = x2;
    i = i + 1;
    % Проверка условия останова
    if sqrt(gr1^2 + gr2^2) <= d;
        break; % Выход из цикла в случае выполнения условия
    end
    k = k + 1;
end
 
Xm = [x1, x2]; % координаты минимума
if k >= kmax
    disp('Минимум не найден!')
    disp('Измените значение шага')
end
 
% Построение графика
d1 = abs(x0(1)-x1);
d2 = abs(x0(2)-x2);
x = linspace(x1-d1,x1+d1,50);
y = linspace(x2-d2,x2+d2,50);
[X, Y] = meshgrid(x, y);
Z = fun(X,Y);
[C, h] = contour(X, Y, Z);
clabel(C, h);
% Отображение меток на линиях уровня
hold on;
plot(x1trace, x2trace, '-');
% Вывод начальной точки на график
text(x1trace(1) + 0.2, x2trace(1) + 0.5, 'M0');
% Вывод решения на график
text(x1 + 0.3*d1, x2, ...
    strvcat(['x1 = ' (num2str(x1))], ...
    ['x2 = ' (num2str(x2))], ...
    ['k = '  (num2str(k))]));
hold off
end

А вот пример использования:

Matlab M

>> f = @(x,y) -exp( -x.^2 - y.^2 - 2*x );
>>GaussSeidelMin(f, [1 1], 1e-3, 0.1)
 
ans =
 
   -0.9976    0.0009

@Сергей.MS 12 / 12 / 1 Регистрация: 16.03.2012 Сообщений: 41
	13.01.2013, 23:42 [ТС]	2
	по читал инфу: итог - нужно в gr1 и gr2 прописать производные соответственно по x1 и x2 от функции f(x1,x2) производные нашёл, прописал: gr1 = c1 + c12x2 + 2c11x1 gr2 = c2 + c12x1 + 2c22x2 но всё равно пока ни как, график должен сходиться в точку, ну а у меня расходится (значения gr1 и gr2 должны уменьшаться) если кто волокёт в этом деле помогите((( 0

Зосима 5221 / 3552 / 372 Регистрация: 02.04.2012 Сообщений: 6,458 Записей в блоге: 17
	18.01.2013, 14:35	4
	Сергей.MS, супер! Особенно график! А как этот метод называется? 0

@Сергей.MS 12 / 12 / 1 Регистрация: 16.03.2012 Сообщений: 41
	18.01.2013, 14:37 [ТС]	5
	Метод покоординатного спуска (Гаусса-Зейделя) 0

Зосима 5221 / 3552 / 372 Регистрация: 02.04.2012 Сообщений: 6,458 Записей в блоге: 17
	18.01.2013, 14:51	6
	А в чем там загвоздка была? Как переписать ее для другой ф-ции? 0

@Сергей.MS 12 / 12 / 1 Регистрация: 16.03.2012 Сообщений: 41
	18.01.2013, 14:56 [ТС]	7
	Вообщем сходил на консультацию в универ, и вот что узнал вся замануха вот в чем: 1. необходимо правильно найти производные от функции двух переменных gr1 b gr2, не сразу получилось) 2. подобрать самостоятельно постоянную шага, так чтобы x1'' было меньше x1' и x2'' меньше x2'. P.S. зря переименовал тему, мне ещё 3 метода решить надо будет)) 1

Зосима 5221 / 3552 / 372 Регистрация: 02.04.2012 Сообщений: 6,458 Записей в блоге: 17
	18.01.2013, 15:05	8
	Спасибо! Сообщение от Сергей.MS мне ещё 3 метода решить надо будет)) Они будут каждый в своей теме! 0

@Сергей.MS 12 / 12 / 1 Регистрация: 16.03.2012 Сообщений: 41
	18.01.2013, 21:59 [ТС]	9
	по ходу возник вопрос, а как изменить код, чтобы график выводился на экран с одинаковым масштабом по осям, если одинаковые величины откладываются, см график. по x -+10 и по y -+10 а в итоге на экране прямоугольник? 0

Зосима 5221 / 3552 / 372 Регистрация: 02.04.2012 Сообщений: 6,458 Записей в блоге: 17
	18.01.2013, 22:43	10
	после plot напиши axes equal 0

@mavi93 0 / 0 / 0 Регистрация: 08.11.2013 Сообщений: 10
	22.05.2014, 00:06	11
	Вот вы сделали с постоянным шагом, а по методу оптимизации Гаусса-Зейделя шаг этот зависит от предыдущих итераций, как это сделать вот ссылочка теорию и пример 0

@olegar4 0 / 0 / 1 Регистрация: 07.06.2013 Сообщений: 29
	26.11.2014, 00:38	12
	этот код работает только для определенной функции, так? А можно создать М-файл для общего решения? 0

@mr_Adrian 0 / 0 / 0 Регистрация: 17.10.2014 Сообщений: 12
	15.12.2014, 01:03	16
	спасибо!. признателен. 0

@droi 0 / 0 / 0 Регистрация: 14.12.2014 Сообщений: 4
	20.12.2014, 20:24	19
	Зосима, Спасибо! 0

@mr_Adrian 0 / 0 / 0 Регистрация: 17.10.2014 Сообщений: 12
	12.12.2014, 02:57	14
	Сообщение от Зосима gr1 = ( fun(x1+h,x2)-fun(x1-h,x2) )/(2*h); % производная по х1 Поскажите , пожалуйста по какому принципу находиться производная.(условно говоря расшивровать эту строчку) 0

Зосима 5221 / 3552 / 372 Регистрация: 02.04.2012 Сообщений: 6,458 Записей в блоге: 17
	12.12.2014, 10:49	15
	mr_Adrian, по определению производная это предел отношения приращения ф-ции к приращению аргумента. Но в данном случае я нахожу не предел, а отношение при очень малом приращении h это простейшая формула численного расчета производной 0

@droi 0 / 0 / 0 Регистрация: 14.12.2014 Сообщений: 4
	20.12.2014, 16:43	17
	Не подскажите стандартные функции оптимизации в matlab? 0

Зосима 5221 / 3552 / 372 Регистрация: 02.04.2012 Сообщений: 6,458 Записей в блоге: 17
	20.12.2014, 17:55	18
	droi, это fminsearch и fminbnd 0

Оптимизация методом покоординатного спуска (Гаусса-Зейделя)

Решение