xset(’window’,1)

plot2d(w,PY(1:lw))

xtitle(’S(w)’)

Результаты приведены на рис. 6.8

Глава 7

Решение задач оптимизации

Оптимизация как выбор наилучшего варианта среди некоторого множества подразумевает наличие правил предпочтения одного варианта другому – критерия оптимальности. В основе построения такого правила лежит целевая функция (функция качества), аргументами которой являются управляемые параметры – внутренние

параметры , которые можно изменять на данном этапе проектирования. Целевую функция обозначим через ( ), область ее опреде-

ления – через . Если область определения дискретная, то задача относится к задаче дискретного (в частном случае целочисленного) программирования. Множество ( ) точек, которые находятся от

точки 0 на расстоянии, не превышающем заданное > 0, называют -окрестностью точки 0 :

( ) = { : ‖ − 0‖ ≤ } .

Максимумом функции ( ) называют ее значение ( *), если суще-

ствует число > 0, такое, что для любой точки ( *) ( ̸= *) выполняется неравенство

( ) − ( *) < 0.

Точку * называют локальным экстремумом (экстремальной точ-

кой). Если функция ( ) имеет один экстремум, то она одноэкстремальна, и многоэкстремальна, если имеет более одного максимума

87

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

< < ,

а область в пространстве параметров, удовлетворяющую прямым ограничениям, называют допустимой областью:

Кроме прямых ограничений есть также функциональные, имеющие вид равенств или неравенств:

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

В качестве примера рассмотрим задачу нахождения экстремума функции двух переменных и ее решения с помощью пакета Scilab.

xref=[1;2;3];x0=[1;-1;1] deff(’[f,g,ind]=cost(x,ind)’,...

’f=0.5*norm(x-xref)^2,g=x-xref’); [f,xopt]=optim(cost,x0) //Simplest call

xopt = 1. 2. 3.

f=

0.

Функция optim имеет следующие основные параметры: [f [,xopt

[,gradopt [,work]]]]=optim(costf ,x0 [,algo] [,df0 [,mem]]

[,work] [,<params>] [,imp=iflag]) Для использования функции optim необходимо определить функцию потерь (целевую функ-

цию), которая имеет следующий вид: [f,g,ind]=costf(x,ind) , где f и g – значения функции потерь и ее градиента, переменная ind используется как входной и выходной аргумент и определяет

настройку вычисления функции потерь. Если ind равен 2, 3 или 4,

то функция должна вычислять только f, g и f или только g. Если ind=0, то функция costf не вычисляет ничего и используется только для вывода.

x0 – вектор вещественных чисел, определяющий начальное значение оптимизируемой переменной;

f – значение функции потерь после оптимизации

(f=costf(xopt))

xopt – найденное лучшее значение x;

algo – строка, с помощью которой задается алгоритм оптимизации (’qn’ – квази-ньютоновский, ’gc’ – сопряженных градиентов или ’nd’ – без дифференцирования);

df0 – приблизительное убывание f после первой итерации (по умолчанию df0=1);

mem – количество переменных, используемых для аппроксимации гессиана, (для algo=’gc’ или ’nd’). По умолчанию – 6;

nap – максимальное количество вызовов фунции costf (по умолчанию – 100).

iter – максимальное количество итераций (по умолчанию –

100).

Для облегчения использования функции optim существует функция NDcost, служащая для вычисления градиента с использованием конечных разностей. Параметры этой функции:

[f,g,ind]=NDcost(x,ind,fun,varargin)

x – вектор или матрица вещественных чисел; ind – целочисленный параметр (см. optim);

fun – функция, имеющая следующий вид: F=fun(x,varargin) , varargin может использоваться для передачи параметров

p1,...pn;

f – значение критерия в точке x (см. optim); g – значение градиента в точке x (см. optim).

Рассмотрим пример минимизации функции Розенброка. Представим ее в системе Scilab:

89

function f=rosenbrock(x,varargin) p=varargin(1)

f=1+sum( p*(x(2:$)-x(1:$-1)^2)^2 + (1-x(2:$))^2) endfunction

Теперь минимизируем данную функцию при помощи функции optim, передавая вместо функции costf список, содержащий обращение к функции NDcost, обращение к пользовательской функции

rosenbrock и параметр varargin:

x0=[1;2;3;4];

[f,xopt,gopt]=optim(list(NDcost,rosenbrock,200),x0)

Получим следующие результат:

gopt = 0.0000002 0.0000002

-1.414D-08

-5.112D-08 xopt =

1.

1.

1.

1.0000000 f =

1.

//Модель линейного колебательного звена:

//x’’(t) + c x’(t) + k x(t) = 0 ,

//

//Перепишем ее в виде системы

//дифференциальных уравнений

//первого порядка, обозначив

//y(1) = x и y(2) = x’:

//

//dy1/dt = y(2)

//dy2/dt = -c*y(2) -k*y(1).

//Пусть есть m измерений of x (соответствующих y(1))

//в разные моменты времени

//t_obs(1), ..., t_obs(m),

//обозначенные x_obs(1), ..., x_obs(m)

//и мы хотим идентифицировать параметры

//c и k минимизируя сумму квадратов ошибок

//между x_obs и y1(t_obs,p).

//

function dy = DEQ(t,y,p)

// Правые части нашей системы ОДУ c =p(1);k=p(2) dy=[y(2);-c*y(2)-k*y(1)]

endfunction

function y=uN(p, t, t0, y0)

//Численное решение получим по ode. y = ode(y0,t0,t,list(DEQ,p))

endfunction

function r = cost_func(p, t_obs, x_obs, t0, y0)

//Это минимизируемая функция,

//являющаяся суммой квадратов ошибок

//между значениями выхода модели

//и результатами измерений

sol = uN(p, t_obs, t0, y0) e = sol(1,:) - x_obs

r = sum(e.*e) endfunction

// Данные

y0 = [10;0]; t0 = 0; // Начальное значение y0

// для начального момента t0. T = 30; // Конечное время для измерений.

// Эдесь мы моделируем данные эксперимента pe = [0.2;3]; // Точные значения параметров m = 80; t_obs = linspace(t0+2,T,m);

// Моменты наблюдения

Least square fit to identify ode parameters