Законы распределения оценок параметров динамических объектов при идентификации методом экспоненциальной модуляции Текст научной статьи по специальности «Математика»

 АНАЛИЗ И СИНТЕЗ СИСТЕМ УПРАВЛЕНИЯ

1 УДК 621.398-501.71

ЗАКОНЫ РАСПРЕДЕЛЕНИЯ ОЦЕНОК ПАРАМЕТРОВ ДИНАМИЧЕСКИХ ОБЪЕКТОВ ПРИ ИДЕНТИФИКАЦИИ МЕТОДОМ ЭКСПОНЕНЦИАЛЬНОЙ МОДУЛЯЦИИ

Д.Н. Анисимов, А.В. Хрипков

Московский энергетический институт (технический университет)

Рассмотрен новый метод идентификации динамических объектов, получивший название метода экспоненциальной модуляции. Показано, что он прост в реализации и обладает высокой помехоустойчивостью. Рассмотрены вопросы, связанные с выбором оптимального входного сигнала. Определены законы распределения оценок при наличии помех типа «белый шум» и «синусоида со случайной фазой».

ВВЕДЕНИЕ

Один из важных этапов проектирования систем автоматического управления заключается в идентификации динамических объектов. Существует немало методов идентификации, например, метод наименьших квадратов в различных модификациях, метод максимального правдоподобия, метод инструментальной переменной, метод модулирующих функций и др. У каждого из них свои достоинства и недостатки. Некоторые методы дают несмещенные оценки параметров объекта, но не эффективны в условиях высоких помех. Другие требуют, чтобы закон распределения помехи имел определенный вид. Иными словами, не существует универсальных методов идентификации, которые полностью удовлетворяли бы разработчика систем автоматического управления.

В статье рассмотрены некоторые аспекты использования метода экспоненциальной модуляции, впервые описанного в работе [1], для идентификации линейных динамических объектов, отличающегося высокой помехоустойчивостью и простотой реализации. Этот метод также не лишен недостатков. В частности, при наличии помех он дает смещенные оценки. Однако смещение, как правило, незначительно и к тому же может быть уменьшено с помощью специальных процедур.

1. МЕТОД ЭКСПОНЕНЦИАЛЬНОЙ МОДУЛЯЦИИ

Пусть дан устойчивый объект с передаточной функцией

параметры ах, ..., ап, Ь0, ..., Ьт которого необходимо определить, Х(б) и У^) — изображения по Лапласу входного и выходного сигналов.

В момент времени 1 = 0 на вход объекта подается произвольный сигнал х(1) и в течение времени наблюдения от 0 до Тн на выходе регистрируется сигнал у(1). Сформируем функцию вида г(0 = ехр(—1/8) такую, чтобы к моменту времени Тн она затухала бы до значений, много меньших 1. Эту функцию будем называть модулирующей, а коэффициент 8 — постоянной времени модулирующей функции (ПВМФ). Перемножим входной и выходной сигналы на эту функцию и вычислим

площади Я. и Я, под образованными таким обрах У

зом кривыми:

Тн Тн

= I х(Ы№1, Яу = | y(1)z(1)d1.

0 0

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

да да да

Ях = |x(1)z(1)d1 - | х(1)г(1^1 * |x(1)z(1)d1 =

= |х(1)ехр(—1/8^1 = Х(8 *),

0

да да да

Яу = І У(1)ї(№1 - | У(1)і(1^1 * | у(1)і(1^1

Ьт s + ...Ьі s + Ь0 у( ч)

W(s) = 1-0 = У-і-У, (1)

П л

ап s + ... а^ + 1

Х( s У

|у(1)ехр(—1/8^1 = У(8 *)

(2)

0

0

н

да

0

0

н

да

0

В этих выражениях Х(9-1) и Д9-1) — значения изображений по Лапласу входного и выходного сигналов, которые были бы получены при подстановке вместо оператора s вещественного числа 9 .

Как следует из выражений (1) и (2):

— Ът9-т + ... + Ъл 9-1 + Ъ0 Щ9 1) = -В---------------1 ----0 = -1. (3)

П~п А-1 1 S

а„ 9 + ... + «19 + 1 х

п

Выражение (3) устанавливает связь между задаваемой ПВМФ 9, вычисляемыми площадями Ях, Я

х У

и неизвестными параметрами объекта а1, ..., ап, Ъ0, ..., Ът. Очевидно, что для вычисления этих параметров необходимо задать г модулирующих функций с постоянными времени 91, ..., 9. и вычислить г пар площадей Ях1, ..., Яхг, Яу1, ..., Яуг, где г = п + т + 1. В этом случае параметры объекта могут быть определены из решения системы линейных уравнений

1 (4)

А = Г 1©,

где А = [ап ... а1 Ьт ... Ь0]т, © = [-8п ... -8^!]т,

Г

1 8, ... 8і

-1 -d18n-

1 8„ ... 8

П - 1

8!

-d181

81

Если идентификация проводится в идеальных условиях (время наблюдения неограниченно, шаг дискретизации бесконечно мал, помехи отсутствуют), погрешность определения параметров равна нулю независимо от вида входного сигнала и выбора ПВМФ. На практике ни одно из перечисленных условий не выполняется, причем основным источником погрешности является помеха п(0, действующая на объект (будем считать, что она приложена аддитивно к выходу объекта). При ее наличии погрешность зависит как от входного сигнала, так и от постоянных времени модулирующих функций. Значения ПВМФ должны задаваться, исходя из следующих соображений. С одной стороны, их нельзя выбирать слишком малыми, поскольку значения вычисляемых площадей в этом случае будут малоинформативными из-за большого удельного веса случайной составляющей. С другой стороны, они не должны быть слишком большими, иначе модулирующие функции не успеют затухнуть к моменту окончания наблюдения процесса, приближенные равенства в выражении (2) будут несправедливы и при вычислении оценок появится регулярная погрешность. Очевидно, должен существовать некоторый диапазон значений ПВМФ, обеспечивающий лучшие в определенном смысле оценки.

Как показано в работе [2], при идентификации методом экспоненциальной модуляции оптимальным входным сигналом является ступенчатое воз-

действие. Такой сигнал обладает рядом преимуществ, а именно:

• отношение «сигнал/шум» максимально среди всех сигналов заданной амплитуды в предположении, что шум имеет нулевое математическое ожидание (при произвольном законе распределения), а объект является фильтром низких частот;

• появляется возможность компенсировать регулярные составляющие погрешностей оценок, и выбор ПВМФ при этом ограничен только снизу;

• значения площадей Sx в выражении (2) априорно известны, что снижает объем вычислений;

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

2. РАСПРЕДЕЛЕНИЕ ОЦЕНОК ПАРАМЕТРОВ ОБЪЕКТА ПРИ НАЛИЧИИ БЕЛОГО ШУМА

Пусть дан объект второго порядка с передаточной функцией

*(0 = -------------

a2S + «^ + 1

и по результатам наблюдения переходного процесса в интервале времени от 0 до Тн необходимо определить параметры а1 и а2. Согласно методу экспоненциальной модуляции, необходимо сформировать модулирующие функции г1 = ехр(—//91) и г2 = ехр(—//92), которые затухали бы к моменту Тн до значений, много меньших 1, умножить на них выходной сигнал и вычислить площади Я1 и Я2.

Тогда оценки «1 и «2 согласно системе уравнений (4) могут быть определены следующим образом:

а =

1

81 - 82

- 8182

а2 = -------------

а2 = -8---1------8-2-

[8? (8і/Яі - 1) - 82 (82/Я2 - 1)],

[8і(8і/Яі -1) - 82(82/^2 - 1)].

Пусть на выходе объекта действует равномерно распределенный в интервале [—А, А] белый шум п(0 с функцией плотности распределения

/Сл(1)) =

1

при п( 1) є [-А; А ],

2А

0 при п(0 £ [-А; А].

Математическое ожидание этого сигнала М(п) = 0, дисперсия £(п) = А /3.

Обозначим через Щ значения площадей, которые были бы вычислены в отсутствие шума,

Г

03

;, а через А5. = 5. — Я. — погреш-

Я =- 2 ..

а2 + а10/ + 0 ности, обусловленные действием шума. Тогда

0

Пусть съем и обработка информации происходят с шагом А?, настолько малым, что погрешностями численного интегрирования можно пренебречь для любого метода интегрирования. Вычислим значения А5., применяя метод трапеций:

(

А5. = А/1

По/2 + X Пкехр(-Ш/0г) к = 1

= А/

V2 + X Пк4 к = 1

где цк = п(кА?), qi = ехр(—А?/0г).

Из этой формулы следует, что математические ожидания величин А5; равны 0, т. е. М(5;.) = Я.. Определим ковариации погрешностей вычисления площадей:

cov(А51, А52) =

= А/2М

г с.

= О(п)А/

0 1 0 2 0 1 + 0 2

Отсюда можно найти дисперсии величин А5;: ДА51) = ^(п)А?01/2, ДА52) = ^(п)А?02/2.

Эти значения представляют собой дисперсии бесконечной суммы независимых случайных величин, однако вклад каждой величины в суммарную дисперсию неодинаков. Поэтому, строго говоря, закон распределения величин А5; отличен от нормального. В качестве меры близости закона к нормальному можно рассмотреть коэффициенты разложения в ряд Маклорена функции 1п[£;(Х)], где g(k) — характеристическая функция случайной величины А5. (для нормально распределенной случайной величины все эти коэффициенты, начиная с третьего, равны нулю [3]):

Лдк

<*Х> = П I

і

■ехр( j ХА/пМп =

к=0 , к 2А4 -Лд

П

к=0

sin ( ХА /Ад ) ХА/Адк

9_ А /

1п[*(Х)] = X 1п

к = 0 2

sin ( X А /Адк )П к

. ХА/Ад .

(ХА /А)2 (ХА /А)4 (ХА/А) 6 2•6 4 • 180 6•2835

Видно, что коэффициенты разложения убывают достаточно быстро, поэтому ограничимся лишь первым членом разложения:

1пк(х)] .......)2

-1 а/0. (ХА! 2.

6 '

Данная характеристическая функция соответствует нормальному закону распределения А5. с нулевым математическим ожиданием и дисперсией

ДА5г) = А? 0гА2/6 = Дг|)А/ 0/2.

Таким образом, можно считать, что при идентификации методом экспоненциальной модуляции влияние равномерно распределенного белого шума эквивалентно влиянию нормально распределенного шума. В этом случае закон распределения вектора S = (51, 52) имеет вид:

/(5„ 52) = + 2с12(51

_ JС11С22 - С12

ехр

2 [Сіі(5і - Я,)2 +

где с,, =

С22

2 П

Я1)( 52 — Я2) + С22(52 — Я2)

О(А52)

О (А 51) О(А 52) - ^2(А 51; А 52)

О ( А 51 )

О(А51)О(А52) - cov2(А51, А52)

cov(А5l, А52)

С12 = — 2 .

О (А 51) О(А 52) - cov2 (А 51; А 52)

После подстановки в эти выражения значений дисперсий и ковариации найдем закон распределения вектора а = (а 1, а 2), воспользовавшись формулой [3]

/(а) = /(5( а ))|/( а )|,

где — якобиан,

д( 51> 52 )

/(а) =

д( а 1, а2)

Обозначим Аа. = а. — а і — отклонения оценок параметров от истинных значений. После подстановок и преобразований получим

2

X

(01 + 02)( 5152)

/(а 1, а 2) = -------1----2----^

2п О (п)А /(0102)

7/2

х ехр

В 2 2 "

—2 (^цАа1 + ^12 А а1Аа2 + ^22Аа2)

(5)

эо

30

о

о

о

о

где В —

(01+ 02)2

2

—1

^11 = — +

8 (—1— 2 )

0

(0102)2(01 + 02)

+ —4.

л5

12

+

2—4

01 (0102)“

—4 8 (—! —2 ) 2

2

й?22 — — +

+ —2. 0о

01 (0! 0 2 ) 3 (0! + 0 2 )

Плотности проекций вектора а на оси О 1 и О 2 определяются, соответственно, следующими выражениями:

/1( а 1) — |/(а!, а2 )^а2,

-о

о

/2( а 2) — | /(а 1, а 2 >^а х.

(6)

(7)

На рис. 1 изображен график плотности проекции вектора а на ось а 1 при истинных значениях

параметров объекта а1 — 3 ^ а2 — 2 c , ПВМФ 01 — 1 c, 02 — 1,1 ^ шаге дискретизации А/ — 0,01 c, соотношении «шум/сигнал» 0,25; на рис. 2 — график плотности проекции вектора а на ось а 2 при тех же условиях. На этих же рисунках изображены гистограммы соответствующих случайных величин, полученные по результатам 50 тыс. экспери-

Рис. 1. График плотности проекции вектора a на ось а построенный в соответствии с выражениями (5) и (6) (кривая 1) и гистограмма, полученная экспериментально (ступенчатая кривая 2)

Рис. 2. График плотности проекции вектора a на ось а 2, построенный в соответствии с выражениями (5) и (7) (кривая 1) и гистограмма, полученная экспериментально (ступенчатая кривая 2)

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

3. РАСПРЕДЕЛЕНИЕ ОЦЕНОК ПАРАМЕТРОВ ОБЪЕКТА ПРИ СИНУСОИДАЛЬНОЙ ПОМЕХЕ СО СЛУЧАЙНОЙ ФАЗОЙ

Для многих промышленных объектов характерно наличие помехи синусоидального вида, обусловленное обычно влиянием электрической сети переменного тока. Амплитуда и частота сигнала помехи, как правило, постоянные, а фаза является случайной величиной с равномерным распределением. Определим плотность распределения оценок для этого случая.

Пусть на вход объекта подается единичное ступенчатое воздействие, а к его выходу аддитивно приложено возмущающее воздействие п(?) = = Азт(ш? + ф), где А — амплитуда сигнала, ш — его частота, ф — случайная фаза с плотностью распределения

/(ф)

1/2 п при - п < ф < п ,

0 в остальных случаях.

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

Ж(5) = огтт • (8)

где а - неизвестный параметр объекта, который требуется оценить. В соответствии с методом экспоненциальной модуляции необходимо сформировать одну модулирующую функцию £(?) = ехр(—?/0), умножить ее на выходной процесс у(?) и вычислить площадь 5 под образованной таким образом кривой.

Точное значение площади Я, которое было бы вычислено в отсутствие помехи при ступенчатом входном воздействии, согласно выражениям (2) и (8), определяется как

2

—

0

а + 0

1

6

о

о

^

Рассмотрим распределение погрешности вычисления площадей А5 = 5 — Я, обусловленной действием помехи. Легко показать, что

ASy = J n(t)exp(-t/e)dt =

о

= J* A sin(«t + ф)єхр(-t/e)dt = о

Ae

J

:sin^ + arctg«0).

1 + (ш0)2

Зная плотность распределения фазы ф и зависимость погрешности А5 от фазы, можно определить плотность распределения величины А5:

/(AS) =

і (%)

- A S2

'U + («О)2

где 1(%) =

0 при % < О,

1 при % > О,

% _ A О

Я

jAS j.

ІІ + («О)

Обозначим Aa = a - a, n = a + О, m =

A ( a + 0 ) . Тогда

eV І + («О)2

/(Aa)

1 ( %)

п| І + A a / n| л/m 2 (n + A a)2 - A a'

. (9)

Найдем область определения этой функции. Как следует из выражения (9):

т2(п + Аа)2 — Аа2 > 0 или — т < А а < т.

А а + п

В зависимости от значения т возможны три решения этого двойного неравенства:

тп тп

І) Aa є

2) Aa є

І + m І

m.

при m < І;

2

при m = І;

З) Aa є I-ад,

mn І - m_

u

mn

_ І + m

при m > І.

Очевидно, второй случай возможен лишь при счетном множестве сочетаний значений А, ш и 0, поэтому рассматривать его не имеет смысла. Третий случай является неблагоприятным, поскольку только правая часть области определения содержит желаемое значение Аа = 0. Что касается первого случая, то его всегда можно достичь правильным подбором значения 0. Чем больше это значение, тем компактнее область определения функции плотности распределения и тем меньше разброс оценок параметра относительно его истинного значения. На рис. 3 изображены график плотности распределения и гистограмма оценки параметра а при истинном значении а = 1 с, ПВМФ 0 = 0,9 с, соотношении «шум/сигнал» 0,15.

ЗАКЛЮЧЕНИЕ

Рассмотрен новый метод идентификации динамических объектов, получивший название метода экспоненциальной модуляции. Показано, что он прост в реализации и обладает высокой помехоустойчивостью. Рассмотрены вопросы, связанные с выбором оптимального входного сигнала. Определены законы распределения оценок при наличии помех типа «белый шум» и «синусоида со случайной фазой». Совпадение полученных зависимостей и экспериментально построенных гистограмм свидетельствует о корректности проведенных расчетов.

о

о

Рис. 3. График плотности распределения (кривая 1) и гистограмма (кривая 2) оценки параметра объекта при воздействии синусоидальной помехи со случайной фазой

ЛИТЕРАТУРА

1. Анисимов Д.Н. Идентификация линейных динамических объектов методом экспоненциальной модуляции // Вестник МЭИ. - 1994. - № 2. - С. 68-72.

2. Анисимов Д.Н. Некоторые аспекты применения метода экспоненциальной модуляции для идентификации динамических объектов // Тр. междунар. конф. «Идентификация систем и задачи управления» SICPRO'03 / Ин-т пробл. упр. - М., 2003. - С. 1602-1618.

3. Пугачев В.С. Теория вероятностей и математическая статистика. - М.: Наука, 1979. - 496 с.

Ш (495) 362-74-07,

е-mail: AnisimovDN@mpei.ru, KhripkovAV@mpei.ru

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

А.С. Манделем. □