Элементы математического моделирования в программных средах MATLAB 5 и Scilab (Андриевский Фрадков)

При этом входной сигнал, создающий максимальное возбуждение системы, будет зависеть не только от времени, но и от состояния системы, т.е. иметь вид обратной связи [129, 108]. Заметим, что для линейных систем задача (3.89) - не что иное, как классическая задана Булгакова о накоплении возмущений [75]. Величина оптимума в задаче (3.89) зависит от 7 квадратично. Поэтому естественно ввести характеристику возбудимости системы ограниченным входным сигналом как величину

Е{7) =

(3-90)

где <2(7) ~ оптимальное значение в задаче (3.89). Лля линейных асимптотически устойчивых систем величина (3.90) не зависит от 7, а для нелинейных - представляет собой функцию от 7, которую естественно назвать степенью возбудимости системы.

Решение задачи оптимального управления в общем случае весьма трудоемко даже при использовании эффективных численных методов. Однако для приближенного решения можно воспользоваться известным результатом (см., например, [117]) о возможности приближенной замены оптимального управления на локально-оптимальное, получаемое максимизацией скорости изменения целевого функционала в силу свободной системы (при и = 0). Точность такой замены тем выше, чем меньше амплитуда входа 7.

Лля вычисления локально-оптимального управления пред-

ставим первую часть уравнения (3.88) в виде
F{x,u)	= /(х) + g(x)u			+ Я(х,и),		(3.91)
где /(х) = F(x, 0),	g(x)	= df/du(x,		0),	а остаточный член
R(x}u) имеет высший порядок			малости		по и.	Скорость из-
менения целевого функционала Qt = y2{t)					равна
Q(t) = 2y(Vh)T		(/(х) + g{x)u + Я(х,				и)).
Пренебрегая величиной R(x,u)y			видим, что локально-оптима-
льное значение входа при малых 7 равно
и(х) =	7sign (л(х) Vh{x)Tg{x)) .					(3.92)

Таким образом, для вычисления степени возбудимости системы при малых 7 достаточно подавать на вход системы сигнал (3.92) и измерять достигаемую амплитуду выхода. Это

141

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

Рис. 3.32. Характеристики возбудимости маятника (д = 0.1, и>$ = 10) и системы Луффинга.

Полученные графики можно использовать для оценки устойчивости замкнутой системы с нелинейностью в обратной связи. Условия устойчивости следуют из теоремы о пассивности [152, 134, 79] и аналогичны (3.86):

K f K 9 < 1.

(3.93)

При этом роль максимума амплитудно-частотной характери-

с т и к и и г р а е т максимальная степень возбудимости

Кр = sup £"(7).	(3.94)
7

Величина Кр конечна для так называемых строго пассивных систем (системы с полной диссипацией).

Отметим, что воздействие вида (3.92) создает в системе аналог резонансного режима: для слабодемпфированных систем малое воздействие вида (3.94) приводит к возбуждению больших колебаний выхода и может сообщить системе значительную энергию. Можно показать, что для механических систем со степенью демпфирования (диссипации) g > 0 воздействие (3.92) выводит систему на уровень энергии не меньший, чем 72 /д2 (см. [108]), т.е. для Е(7) справедлива нижняя оценка £ ( 7 ) > д~~1.

142

Исследование динамических свойств систем с помощью непериодических тестовых сигналов представляется весьма перспективным инструментом в теории нелинейных систем.

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

(теория связи). Однако и	здесь	"нелинейная	философия"
предлагает новые подходы.
Напомним, что в теории связи гармонический			сигнал
y(t) =	asin(utf	+ a)	(3.95)

рассматривается как базовый, простейший. Сигнал может изменяться (модулироваться) путем изменения его параметров - величин а (амплитуда), и (частота) и а (начальная фаза). Например, при частотной модуляции модулирующим параметром является частота а;, которая становится, таким образом, переменной, т.е. с системной точки зрения становится сигналом: и = и>(t). Передаваемый модулированный сигнал y(t) может содержать в себе закодированное сообщение. Для передачи сообщения по каналу связи, на стороне приемника сообщение должно восстанавливаться по принимаемому сигналу y(t) = y(t) + £(t), где £(t) - шум (помеха) в канале связи. Выделение полезного сообщения из принимаемого сигнала и является основной задачей теории связи [25].

В 90-х годах XX века усилился интерес к использованию в качестве несущих сигналов нерегулярных (хаотических) колебаний [33, 143, 133, 6]. Для построения соответствующей теории необходимо сделать решительный шаг: перейти от явного описания сигнала как функции времени (3.94) к заданию модели системы, генерирующей этот сигнал. Очевидно, например, что генератором гармонического сигнала (3.94) может служить линейное дифференциальное уравнение

y(t)+u2y(t)

= 0.

(3.96)

Однако для генерации нерегулярных, непериодических сигналов линейные модели непригодны (см. п. 3.8.1). Следующим шагом к построению новой теории является использование в качестве основных объектов нелинейных генераторов,

143

описываемых дифференциальными уравнениями

x =	F{xtu)t	у = h(x).	(3.97)
Здесь в отличие от	(3.88)	входной вектор и £ 1Zm	может

иметь смысл не только входного сигнала, но и модулирую-

щего	сигнала, т.е.	задавать	набор	изменяемых		параметров
генератора. Например, система (3.97)					в частном	случае мо-
жет	быть системой	Лоренца	(3.63)	или	системой	Чуа (3.66),

а вектор входов может включать часть коэффициентов соответствующего уравнения.

Приемник также представляется нелинейной динамической

системой
z = 9(z,y), « = *(*),	(3.98)

где у = 2/(0+ £ (0 - принимаемый сигнал (вход приемника), z = z(t) € TV1 - вектор состояния приемника, й = u(t) - оценка передаваемого сообщения (выход приемника). Обычно шум измерений можно считать ограниченным. Тогда задачу синтеза (конструирования) приемника можно поставить как нахождение модели (3.98), обеспечивающей достижение цели

|u(t) - 4 ( 0 1 < Дм.			(3.99)
На самом деле, соотношение (3.99) может не		выполнять-
ся на начальном этапе работы устройства. Поэтому			цель в
задаче синтеза формулируется как асимптотическое			соотно-
шение
lim	М О " f i ( 0 l <		(3.100)
В частном случае Ди	= 0 цель (3.100) означает	асимптоти-

чески точное оценивание (о нем имеет смысл говорить, если пренебрегать помехами).

Естественным подходом к решению задачи синтеза является включение вектора оценок состояния передатчика в вектор состояния приемника z(t). Например, выбирают z(t) = x(t) £ 6 TZn . Тогда система (3.98) представляет собой не что иное, как наблюдатель (или фильтр состояния) для системы (3.97) (см. [6]).

144

Вэтом случае естественно поставить дополнительную цель

-достижение заданной точности оценивания (наблюдения) состояния передатчика

Tim \x(t) - х ( 0 | < Д..

(3.101)

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

Синтезу адаптивных наблюдателей (фильтров состояния) для линейных систем посвящена обширная литература (см. [6]). При наличии помех их обычно считают случайными и ставят задачу оптимальной фильтрации - достижение минимальной среднеквадратической ошибки оценивания параметров. Однако включение вектора оценок неизвестных параметров в общий вектор состояния системы делает задачу фильтрации нелинейной даже при линейных моделях передатчика и приемника, что приводит к значительным вычислительным трудностям.

Синтез наблюдателей для нелинейных систем в общем случае также представляет собой трудную задачу. Однако для частного класса так называемых "пассивируемых" систем [63, 79] решение оказывается простым и изложено ниже, в п. 6.2. Заметим, что такие генераторы хаоса, как системы Лоренца и Чуа (см. п. 3.8.1), в большинстве случаев являются пассивируемыми.

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

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

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

Карл Поппер

ГЛАВА 4. ВЫБОР П А Р А М Е Т Р О В МАТЕМАТИЧЕСКОЙ М О Д Е Л И

4.1.Предварительные преобразования

4.1.1.Линейно-параметризованные модели

Итак, мы выбрали структуру ММ системы, т.е. выбрали ММ с точностью до конечного набора числовых параметров. Как говорят, "модель параметризована". Что же дальше? Значения параметров лишь в редких случаях удается подобрать исходя только из теории или априорных соображений. Как правило, оценка параметров ММ проводится по результатам наблюдений за реальным процессом или явлением в ходе нормального функционирования либо во время специальных экспериментов. Более того, без сопоставления результатов наблюдений за реальной системой и за ее ММ нет гарантий правильного выбора структуры ММ. Поэтому на этапе выбора параметров происходит окончательное уточнение и окончательный выбор структуры ММ, желательно из нескольких конкурирующих вариантов.

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

N

(4.1)

i=i

где х, - входные переменные (входы, факторы), у - выходная переменная (выход, отклик), б, - параметры ММ, <р- возмущение (погрешность ММ). Если ввести векторные обозначения

146

X = c o l ( x 1 ? х п ) е	Rn}e = col(0b ...,0n ) е Rn	, то (4.1) можно
переписать в виде	у = вТх + <р.	(4.2)

По виду (4.1) ясно, что эта ММ приспособлена для отражения количественных зависимостей, т.е. множества значений входов и выходов системы U, Y должны быть непрерывными.

Этот случай (при случайных <р это задачи регрессионного анализа) будет подробно рассмотрен ниже. Иногда моделью (4.1) можно описать и зависимость между величинами с дискретным множеством значений (например, xi}y могут быть номерами уровней входных переменных, т.е. принимать значения 1,2,3,...). Однако для исследования качественных зависимостей более приспособлены методы дисперсионного анализа [46, 110], которые здесь рассматриваться не будут. Подчеркнем, что величины х,у,0,<р, вообще говоря, не совпадают с исходными входами, выходами, параметрами и возмущениями в системе, а получаются в результате преобразования (замены переменных), специально сделанного для приведения ММ системы к форме (4.1). Покажем, как можно перейти к форме (4.1) для основных типов моделей: статических; динамических, дискретных по времени; динамических, непрерывных по времени.

Лля простоты изложения возмущение в (4.1) учитывать не будем.

4.1.2.Преобразование статических моделей

Линейная М М . Чтобы привести линейную ММ с m вхо-

дами	у = а0 + d\U\ + ... + amum	(4.3)
к виду (4.1), полагаем хх = 1, вх = а0 и для i = 2, ...,га+1, х, =		х,
0, = а,_!.	Таким образом, х = со1(иь ...,um), в = (а0)ах,...	,аш ),
n = rn +	1.

Нормированная линейная М М . Рекомендуется перед оцениванием параметров ММ проводить преобразование нормировки (масштабирование и центрирование) переменных, приводящее диапазоны изменения переменных к стандартно-

му значению.	Это преобразование для линейной ММ соот-
ветствует выбору в (4.1) х,- = — — L .		Обычно выбирают чи-
сла Ui, Gi так,	чтобы нормированные	переменные лежали в

147

диапазоне [—1,1], т.е. щ - "номинальное" или среднее значение переменной, а сг, - максимальное отклонение от среднего.

Встохастическом случае берут щ = Мч{) <7, = у/М(и{ — й,)2. Полиномиальная М М . Полиномиальная ММ с одним

входом и выходом, имеющая вид

у = а0 + ахи + а2и2 + ... + атит}

(4.4)

где т -предполагаемая степень полиномиальной зависимости, приводится к (4.1) выбором

х = col(l,ti, ...,um), в = со1(а0,а!, ...,ат ), п = т + 1.

Квадратичная многофакторная ММ. Квадратичная зависимость с несколькими входами, имеющая вид

оо	оо
у = а0 +	+ ^ dijUiUj,	(4.5)
»=i	»j=i

может быть приведена к форме (4.1), если положить:

в = c o l ( a 0 , a m , a u , ( a 1 2 + a 2 i ) , ( a i m + a m l ) , a22, (a23 +a32 ), •••> ^mm)- Достаточно брать симметричную матрицу старших коэф-

фициентов:

at;

=

a; i .

Общее число

входов в (4.1)

будет

71 =

.

,

,

т ( т

+

1)

=

(ш + 1)(тп +

2)

1 + т Н ^

л

^

L .

Мультипликативная (степенная) ММ . Бели ММ запи-

сана в мультипликативном

виде

у =

Си?

• и?

(4.6)

то для приведения к виду (4.1) следует

прологарифмировать

(4.6):

In у =

In С

+ Qi In Ui + ... +

a m In ит)

и положить у = \пу,

= 1, 01 =1пС и для i = 2 , . . . , г а + 1 ,

= 0ft-«lf

x , = ln 11,-1.

Неявная MM. Если статическая характеристика системы

представлена в виде неявной зависимости

ф ( « , » ^ ) = 0 ,

(4.7)

148

то нет необходимости разрешать (4.6) относительно выхода. Достаточно свести (4.7) к виду

		п
		i=i	= 0,			(4.8)
линейному относительно параметров, и положить
	х = col($i(u,y),..., Фп(<а, у),		у =	Фо(и}у).
	Приведем несколько примеров, демонстрирующих приемы
построения линейно-параметризованных ММ.
	Пример 4.1.1. Одним из распространенных приемов явля-
ется логарифмирование уравнений.			Если ММ задана			экс-
поненциальной зависимостью у = C\eC2U,				то, логарифмируя,
получим соотношение In у = InCi + С2и,			для приведения		кото-
рого к виду (4.1) достаточно положить у = In у, хх = 1, х2						=	и,
в\ = InCi, 02 = С2. После построения			оценок параметров				0 Ь
02	следует вернуться к параметрам исходной ММ: С\				=		е6х,
С2	— $2'
	Пример 4.1.2.	Уравнение Аррениуса,		определяющее			ки-
нетику химических		реакций, имеет вид
	*«*.«ф(-д(Г|273))'				(4,9)
где К - константа скорости реакции,			Т -	температура	в °С,

Е - энергия активации, R - универсальная газовая постоянная. Величина R известна из теории, а Ко, Е подлежат определению из опытов. Прологарифмированное уравнение (4.9)

сводится к (4.1)	подстановкой:
у = In К, xi = 1, х2 =				>	~ 1п ^о»	^2 =	Е.
Пример 4.1.3. Тригонометрическая					модель
	у =	y4sin(u>* + <р),					(4.10)
в которой си, t	известны,	Л,	(р неизвестны, сводится				к (4.1),
если записать	(4.10) в виде		у	= Л cos (р sin u>t +		Л sin <р cos ut и
обозначить: хх	= sin ut, х2 =			cos ut,	= Acosip,	в2 =	Л sin у?.

149

Обратный переход после построения оценок			в2 выполняет-
ся по формулам
	А = в2х+в1	ф = arctgT2-.
П р и м е р	4.1.4. Пусть в предыдущем примере и			неизвест-
на. Тогда	следует сначала	оценить и>, исключив		из (4.10)

остальные параметры. Это можно сделать, приведя двукратным дифференцированием (4.10) к виду у + ш2у = 0 и оценив и2 одним из методов, описанных ниже в п. 3.1.3. После этого можно применить подход предыдущего примера.

П р и м е р 4.1.5. Гиперболическая зависимость у =			^ + fiu
преобразуется	к виду	^ = ^ + /3, который сводится	к (4.1)
заменой у =	хх =	вх = /3, в2 = а.

4.1.3.Преобразование д и н а м и ч е с к и х моделей

Проще всего обстоит дело с дискретными по времени				линей-
ными моделями, заданными передаточной функцией W(z~l)					=
= B(z~l)/A(z~l),	где A(z~l)	= 1 + a1 z"1 + ... + атуг—«<, B(z'1)			=
= bQ + bxz+	... -I- bm%z~m*}	или соответствующим разностным
уравнением
Ук + а\Ук-\ +	••• + а>туУк-ту	= Ь0ик	biUk-x + ... + ЬтчЩb-mu-	(4.11)

ММ (4.11) непосредственно приводится к (4.1), если положить

п =		ти + ту +1;		у = ук,
X	=	Со1( Ук — \)		Ук—тпу I	,Ujt_mJ,
в	=	col(ab	...,аШу,	,Ь0,. -,Ьши).
Если ММ задана уравнениями				состояния
		хк+х =	Ахк + Ьик) у = Схк)		(4.12)

то следует привести их к виду (4.11), как это описано в п. 2.2.2. Уравнение (4.11) содержит меньше неизвестных параметров и включает только входные и выходные, т.е. измеряемые, переменные. Если же измерению доступен непосредственно вектор состояния хк 6 Я", то целесообразно представить векторное уравнение (4.12) как совокупность из п покомпонентных скалярных уравнений, коэффициенты которых

150