ISSN 1814-1196 Научный вестник НГТУ том 80, № 4, 2020, с. 47-64
http://journals.nstu.ru/vestnik Science Bulletin of the NSTU Vol. 80, No. 4, 2020, pp. 47-64
ИНФОРМАТИКА, ВЫЧИСЛИТЕЛЬНАЯ ТЕХНИКА И УПРАВЛЕНИЕ
INFORMATICS, COMPPUTER ENGINEERING AND CONTROL
УДК 681.51:519.6 Б01: 10.17212/1814-1196-2020-4-47-64
Алгоритмы непараметрической идентификации
*
сложных технических систем
Ю.Е. ВОСКОБОЙНИКОВ1'2'", В.А. БОЕВА1'4
1 630008, РФ, г. Новосибирск, ул. Ленинградская, 113, Новосибирский государственный архитектурно-строительный университет 630073, РФ, г. Новосибирск, пр. К. Маркса, 20, Новосибирский государственный
технический университет
а voscob@mail.ru ь v.boyeva@sibstrin.ru
На практике часто встречаются сложные технические системы, представляющие собой соединение нескольких разнотипных более простых подсистем. Из-за сложности физических процессов, протекающих в них, невозможно построить подробную (детальную) математическую модель, адекватно описывающую процесс в каждой подсистеме. В таких случаях берется модель «черного ящика», «внутренности» которого не детализируются. Для стационарных линейных систем (подсистем) в качестве соотношения, устанавливающего связь между входом и выходом «черного ящика», принимается интегральное уравнение Вольтера первого рода с неизвестным разностным ядром, которое в теории автоматического регулирования называется импульсной переходной функцией системы (ИПФ). Поэтому для использования модели «черного ящика» необходимо оценить эту неизвестную ИПФ. Это задача непараметрической идентификации, и для сложных систем ее необходимо решить как для всей системы в целом, так и для каждой подсистемы в отдельности, что существенно усложняет процедуру идентификации. Формально оценивание ИПФ можно рассматривать как решение интегрального уравнения первого рода относительно его ядра по зарегистрированным (с погрешностями) дискретным значениям входного и выходного сигналов. Такая задача является некорректно поставленной, поскольку решение может обладать неустойчивостью относительно погрешностей (шумов измерения) исходных данных. Для получения единственного и устойчивого решения используют регуляризирующие алгоритмы, но специфика входных и выходных сигналов в эксперименте по идентификации ИПФ не позволяет использовать их вычислительные методы (СЛАУ или дискретное преобразование Фурье). Поэтому в данной работе для решения задачи идентификации сложных систем предлагаются два алгоритма идентификации, которые в полной мере учитывают специфику решаемой задачи. В этих алгоритмах оценки ИПФ строятся с использованием первых производных от сигналов идентифицируемой системы, для устойчивого вычисления которых применяется сглаживающий кубический сплайн с выбором параметра сгла-
Статья получена 14 августа 2020 г. Исследование выполнено при финансовой поддержке РФФИ в рамках научного проекта № 20-38-90041. Acknowledgments: The reported study was funded by RFBR, project number 20-3890041.
живания. Приводятся результаты идентификации сложной системы «воздухонагреватель-вентилятор-помещение», показавшие эффективность предлагаемых алгоритмов.
Ключевые слова: задача непараметрической идентификации, интегральное уравнение Вольтерра I рода, интегральное уравнение Вольтера II рода, некорректно поставленные задачи, сглаживающие кубические сплайны, устойчивые алгоритмы идентификации, оценивание оптимального параметра сглаживания, решение практической задачи идентификации сложной системы
1. ВВЕДЕНИЕ И ПОСТАНОВКА ЗАДАЧИ
На практике часто встречаются сложные технические системы, представляющие собой соединение (последовательное или параллельное) нескольких разнотипных более простых подсистем с разными физическими процессами в них. Таким примером может служить теплофизическая система «воздухонагреватель-вентилятор-помещение», основная задача которой -стабилизация параметров микроклимата в некотором помещении (объект регулирования). Нестабильность этих параметров возникает вследствие изменяющихся тепловых возмущений в помещении и воздухонагревателе, реакция системы на которые определяет качество параметров микроклимата и устойчивость системы регулирования [1, 2]. На рис. 1 изображена блок-схема этой системы.
Рис. 1. Блок-схема системы регулирования микроклимата
Fig. 1. A block-diagram of an environmental system
Из-за сложности физических процессов, протекающих в каждой из трех подсистем, невозможно построить детальную математическую модель, адекватно описывающую процесс в каждой подсистеме. Поэтому была взята модель «черного ящика» с входным и выходным сигналами. Был проведен обширный физический эксперимент с разными режимами работы, подтвердивший гипотезу о линейности и стационарности как каждой подсистемы, так и системы «воздухонагреватель-вентилятор-помещение» в целом [3]. Поэтому в качестве характеристики, достаточно полно описывающей динамику, была принята импульсная переходная функция (ИПФ), а в качестве математической модели в терминах «вход-выход» взято интегральное уравнение Вольтерра первого рода с разностным ядром [4]. Для рассматриваемой системы регулирования в целом (см. рис. 1) уравнение (т. е. математическая модель) имеет вид
(t -х) 9(x)dт = f (t), t е[0,Г],
(1)
где к (1) - импульсная переходная функция (ИПФ) системы (ядро интегрального уравнения); ф(т), /(1) - входной и выходной сигналы всей системы в целом. При этом выполняется условие к (1) = 0 при 1 < 0, которое определяется технической реализуемостью системы. Аналогичные интегральные уравнения были приняты для описания динамики каждой подсистемы:
\к1 (1 -т)фу (т)ёт = /у (1), t е[0,Г ], (2)
0
где у = 1, 2, 3 - номера подсистем. Таким образом, для моделирования и управления системами с математическими моделями (1), (2) необходимо оценить (вычислить) ИПФ как всей системы в целом, так и каждой подсистемы. Такая задача получила название задачи непараметрической идентификации, и она сводится, как правило, к решению интегрального уравнения (1) (или (2)) по зарегистрированным (с определенными погрешностями) значениям входного и выходного сигналов идентифицируемой системы. На рис. 2 приведена схема проведенного эксперимента по идентификации рассматриваемой системы, в котором на вход системы подается в определенный момент времени 11 (принимается 11 = 0) прямоугольный импульс вида
[0, если 1 < ?1,
ф н (1) = \ (3)
[ А, если 1 >
Если амплитуда A = 1 (это будет подразумеваться в дальнейшем), то функция (3) называется функцией Хэвисайда. Выходной сигнал fiH (t) подсистемы 1 (на рис. 2 обозначено ПС1) является входом подсистемы 2 (ПС2) Ф2^), а выход fi(t) ПС2 - входом подсистемы 3 (ПС3) фз^). Выходные сигналы, являющиеся реакцией на входную функцию Хэвисайда, имеют нижний индекс Н - это функции fi^ (t), fn (t).
Рис. 2. Схема эксперимента по идентификации ИПФ
Fig. 2. An experiment of impulse response identification
Из-за погрешностей измерений в эксперименте регистрируются в моменты времени (необязательно равноотстоящие) , 7 = 1...N, «зашумлен-ные» значения:
Ли, = Ли (tí) + ),
~/2г = /2(4 ) +Л2(',"), (4)
и = /з, = Ык) + ЛзС,"), , = 1. N,
где погрешности измерений ,), ^2(^7), ) являются случайными ве-
2 2 2
личинами с нулевыми средними и дисперсиями с^ , с , с соответственно. Эти зарегистрированные значения (4) являются исходными данными для идентификации ИПФ к (т), kj (т) (см. выражения (1), (2)).
Заметим, что задача непараметрической идентификации относится к некорректно поставленным задачам (НПЗ), где решение может не существовать, быть не единственным и не устойчивым к погрешностям задания исходных данных [5].
Для вычисления устойчивого единственного решения НПЗ используются различные (детерминированные или статистические) регуляризирую-щие алгоритмы (РА) [5, 6]. Для применения РА исходное уравнение (1) или аппроксимируют системой линейных алгебраических уравнений с плохо обусловленной матрицей, или заменяют дискретной сверткой и строят РА на основе дискретного преобразования Фурье (подробнее см. [6, 7]). К сожалению, эти подходы к построению РА неприемлемы для рассматриваемой задачи идентификации. Во-первых, два входных сигнала (для ПС1 и всей системы в целом) являются функцией Хэвисайда, и это вызывает проблемы при применении дискретного преобразования Фурье. Во-вторых, входные сигналы для ПС2, ПС3 содержат случайные погрешности, которые эффективно не учитываются при построении РА вышеназванными подходами (в лучшем случае при выборе параметра регуляризации - принцип обобщенной невязки).
Поэтому основной целью данной работы является построение устойчивых алгоритмов непараметрической идентификации сложных систем с различными типами экспериментальных сигналов. Для этого решаются две задачи.
Задача 1. Построение эффективного алгоритма непараметрической идентификации, когда входной сигнал идентифицируемой системы задается функцией Хэвисайда (обозначим его как алгоритм идентификации 1, или АИ1).
Задача 2. Построение эффективного алгоритма непараметрической идентификации для произвольного входного сигнала идентифицируемой системы (обозначим как АИ2).
При решении этих двух задач предполагается, что соответствующие сигналы регистрируются со случайными погрешностями - шумами измере-
ний (см. (4)). Это обусловливает при вычислении первых производных сигналов (используемых в алгоритмах идентификации) применение сглаживающих кубических сплайнов как средства для эффективной фильтрации шумов измерений. Основные понятия, необходимые для применения сглаживающих кубических сплайнов в алгоритмах идентификации, будут приведены позже.
2. АЛГОРИТМ НЕПАРАМЕТРИЧЕСКОИ ИДЕНТИФИКАЦИИ 1 (АИ1)
Предположим, что идентифицируемая система является стационарной и в качестве математической модели принято интегральное уравнение (1). Во многих практических схемах идентификации ИПФ такой системы на ее вход подается функция Хэвисайда. Для такого входного сигнала можно показать [4], что
к(0 = ~/и (0, ( е[0,Т], (5)
&
где /и (() - выходной сигнал (реакция системы), если на вход подана функция Хэвисайда. Несмотря на хорошо разработанные алгоритмы численного дифференцирования, использование (5) на практике связано с той же проблемой некорректности, что и при решении уравнения (1), так как операция дифференцирования является некорректно поставленной задачей [5]. Чаще всего это проявляется в неустойчивости операции дифференцирования, когда даже небольшой уровень шума регистрации выходного сигнала вызывает очень большие ошибки в полученной (на основе (5)) оценке ИПФ. Для получения устойчивой оценки ИПФ зарегистрированный (с ошибками измерения)
сигнал /и (^) необходимо первоначально сгладить (т. е. отфильтровать шумы), а затем применить операцию дифференцирования. Часто на практике для этих целей используют (из-за их простоты построения и наличия соответствующего программного обеспечения в ряде математических пакетов) сглаживающие кубические сплайны (СКС) [8-10].
Предположим, что зарегистрированный в узлах 0 = ¿1 < ^ <. • • < tN = Т
зашумленный выходной сигнал /и допускает представление (см. (4)): /и! = /и (¿7) + Ц, 7 = 1.N, где ц - шум с нулевым средним и дисперсией сЦ. Для сглаживания зашумленных значений |/и. | обратимся к СКС.
Напомним, что СКС S/a (^) на каждом отрезке , ¿7+1) , 7 = 1. N -1, представляет собой полином третьей степени вида [8, 9]
я/,«(х)=а + Ь С - Ч)+« - ь )2 + &{« - )3 (6)
и является дважды непрерывно дифференцируемым на всем интервале [0,Т]. Для однозначного вычисления коэффициентов СКС а{, Ъ^, С7, &
задают левые и правые краевые условия [9, 10]. Было показано [11], что в задаче параметрической идентификации целесообразно задавать не естественные краевые условия вида Бу,а (/1) = Бу-,а(¿и) = 0 [12], а первые производные:
Б'у^Ь) = .1, Уа(^) = N, (7)
где .1, .И - задаваемые значения (чаще всего нулевые) идентифицируемой ИПФ на концах интервала [0,Г]. Показано [8, 13], что СКС с такими краевыми условиями доставляет минимальное значение функционалу:
N |2 п , _ 0
^ (Б) = а | у (0| Л + Х р,- (/и, - Бу,а (Ч))2 , (8)
/1 ,=1
где р- \ - = 1... N - весовые множители (задаются пропорционально дисперсии шума измерения с ). Как видно из (8), СКС (в отличие от интерполяци-
Лг
онного) не проходит через значения /и1, и поэтому он используется для
фильтрации (сглаживания) зашумленных значений. Параметр сглаживания а «управляет» гладкостью сплайна (а следовательно, и ошибкой сглаживания), и, меняя параметр сглаживания а в интервале [0, да), можно варьировать ошибку сглаживания. Так, при а = 0 сглаживающий сплайн становится интерполяционным (т. е. проходит через все заданные точки /и., и фильтрация
шума измерения отсутствует). При а ^да СКС вырождается в прямую линию (эффект переглаживания зашумленных данных). Величину параметра сглаживания, минимизирующего среднеквадратическую ошибку (СКО) сглаживания [13], назовем оптимальным параметром и обозначим а0рГ. Позже будет рассмотрен один алгоритм оценивания а0рГ на основе проверки
статистических гипотез об оптимальности того или иного значения параметра сглаживания. Вычисленная оценка (обозначим ее как а^) с приемлемой точностью оценивает а0рГ, и в дальнейшем предполагается, что СКС строится при а = а^ .
Сам алгоритм вычисления коэффициентов СКС при заданном а и краевых условиях (7) подробно изложен в работах [8, 13] и здесь не приводится. Только заметим, что первоначально из решения системы линейных алгебраических уравнений с пятидиагональной матрицей находится вектор значений второй производной СКС в узлах сплайна, через который затем вычисляются все коэффициенты, входящие в представление (6).
Таким образом, нахождение оценки ка (т) на основе соотношения (5)
сводится к построению СКС по зарегистрированным значениям /и.,
1 = 1... N, с краевыми условиями (7) и вычислению первой производной СКС, т. е.
4 со=—^ со. (9)
— х
Из (6) непосредственно следует выражение для первой производной СКС: — 2
—3/,а№ (х) = 3'^ (х) = Ъ + 2су (х - Ч) + 3—1 (х - Ч )2, (10)
—х
если ^ <х < 7^+1.
В качестве иллюстрации эффективности работы изложенного алгоритма идентификации приведем некоторые результаты решения задачи идентификации ПС1 (подсистема «Воздухонагреватель»). На вход этой системы подавалась функция Хэвисайда (3), зарегистрированные значения Лн , 1 = 1...N,
показаны на рис. 3 точечной кривой. Следует заметить, что выходной сигнал Лн (7) имеет запаздывание (по сравнению с входным сигналом), равное х г1 = 18 мс, но это запаздывание на рис. 3 не показано.
Очевидно, что наличие значительных погрешностей (относительный уровень примерно 5 %) приведет к большим ошибкам идентификации при вычислении производной по интерполяционному сплайну (который не осуществляет фильтрацию шума измерения выходного сигнала). На рис. 4 точечной кривой показаны значения такой оценки ко (х). Видны осцилляции, обусловленные шумами измерения выходного сигнала и являющиеся характерной приметой неустойчивого решения некорректно поставленной задачи идентификации ИПФ.
Для получения устойчивой оценки ИПФ по этим же исходным данным была разработана следующая методика. Если значение первой производной справа можно обоснованно задать равным нулю, т. е. 3'/ (^) = 0 , то значение к (0) априори не известно. Поэтому первоначально строился СКС со смешанными краевыми условиями 3/ (0) = 0, 3/(tN) = 0 и параметром сглаживания а = а^, первая производная этого сплайна бралась в качестве «стартовой» оценки к^ (х). Затем строился набор из десяти СКС, у которых менялись только значение левого краевого условия (из некоторой
окрестности значения к0° (0)), и из этого набора «регуляризированных»
а^
ИПФ выбиралась «наилучшая» оценка к (х), имеющая наименьшую вели-
а^
N , „ _ ч2 „ „
чину функционала ДА« ) - Лш. ] , где с Иа№ (71) = |ка№ (х) — х - про-1=1 0 гнозные значения выходного сигнала ПС 1.
Рис. 3. Выходные сигналы ПС1 Fig. 3. Output signals of subsystem 1
На рис. 4 сплошной кривой показаны значения такой оценки каw (т).
Видно отсутствие шумовых осцилляций и хорошее согласие с априорным представлением об идентифицируемой ИПФ как о достаточно гладкой функции. На рис. 3 сплошной кривой показаны значения (ti), вычисленные
по наилучшей к (т). Видно, что такой «прогнозный» выходной сигнал
представляет собой «хорошо» сглаженные значения измеренного выходного сигнала, что косвенно подтверждает адекватность построенной ИПФ экспериментальным данным.
0.04 0.02 0 0.02 -0.04
'л-К : ' ...... ' 'ill I V./.-IM
„Ч '-■■.''Vi4
. ™ ¡^'.'Л'^1.—
•х. - ',' Л »w.«
t, мс
100
200
300
Рис. 4. Оценки ИПФ ПС1
Fig. 4. Impulse response estimates in subsystem 1
Этот факт и характер вычисленной оценки ИПФ ка^ (т) позволяет говорить об эффективности предложенного алгоритма идентификации ИПФ, когда на вход идентифицируемой системы подается прямоугольный сигнал (типа функции Хэвисайда).
3. АЛГОРИТМ НЕПАРАМЕТРИЧЕСКОЙ ИДЕНТИФИКАЦИИ 2 (АИ2)
Изложим этот алгоритм идентификации, следуя результатам работы [6, 14] и обобщая их для решения задачи идентификации сложных систем.
Дифференцируя уравнение (1) по переменной 7 и выполнив несложные преобразования, приходим к интегральному уравнению Вольтерра II рода
1 г / '(7)
£(0 +-Гф'(Г-т) к(т)^т = ^^-, 7 е[0,Л, (11)
Ф(0) о ф(0)
решение которого уже является корректно поставленной задачей. Однако построение алгоритма оценивания к (7) из этого уравнения сталкивается со следующими основными трудностями: дифференцирование зашумленных входного и выходного сигналов идентифицируемой системы; вычисление интеграла свертки с наименьшими ошибками интегрирования для уменьшения общей систематической ошибки алгоритма идентификации. Для преодоления этих трудностей вновь обратимся к математическому аппарату СКС.
Будем считать, что по зашумленным значениям входного и выходного сигналов идентифицируемой системы
ф,- =ф^) + , / = ¡(и) + л,-, г = 1...N, (12)
где , - случайные погрешности регистрации сигналов с дисперсиями
2 2
о^, о^ соответственно, построены СКС (7), а^ (7) с производны-
ми £ф,а^ (7), 3/ацг (7). Далее обратимся к интегралу свертки, входящему в это уравнение. Для его дискретизации предположим, что на каждом интервале , 7] +1) ИПФ к (7) постоянна и равна к]- = к (/у ), ] = 1. N -1. Тогда интеграл свертки в (11) будем аппроксимировать суммой:
Г
71=0
г-1 ~
ф^ -т) к(т)^ т = Е к]
у =1
V+1
3ф,аж -т)а т
, г = 2.N. (13)
7
Введем в рассмотрение треугольную матрицу Ф' размером (Ж -1) х (Ж -1) с элементами:
Ф' . =
7,З
+1
| ^ф,а ('г+1 -х)^х если ] < ^ <] (14)
0, если ] > г,
где г = 1... N -1, ] = 1... N -1. Штрих в обозначении матрицы указывает, что ее элементы вычисляются через первую производную СКС. Тогда правую часть в формуле (13) можно представить в виде
'г+1 N-1
I (^+1 -т)£(х)аТ=£Ф.]., г = 1,...,N-1, (15)
/1=0 ]=1
где Ф/ ] - г-, .-элемент матрицы Ф'. В дальнейшем предполагается, что шаг дискретизации по времени постоянный и равен Дг, а моменты регистрации = (г - 1)Дг, г = 1,...,N. Тогда можно доказать следующее равенство:
+1 'к+1
('г+1 -х)^ х=| Яф^ (х)^ х, (16)
'] 'к
где к = г - ] +1. С учетом этого равенства элементы матрицы (14) можно переписать в виде
Ф' . =
7,З
'г-]+2
I ^ф а (х)^х, если ] < г, 'г-]+1 (17)
0, если ] > г,
где г = 1... N -1, ] = 1... N -1. Используя тождество (16), интеграл от производной (х) СКС при вычислении элементов матрицы Ф' (см. (17)) можно определить по формуле
'к+1
I ^ (х¥х = ЪкД1 +Ск дД +йкд], (18)
'к
что существенно упрощает вычисление элементов матрицы Ф'.
Заметим, что полученная квадратурная формула (18) позволяет достаточно точно и эффективно вычислять значения интегралов, тем самым уменьшая методическую ошибку предлагаемого алгоритма дискретизации уравнения (11).
Сформируем векторы:
/ ' =
, ка -
"а N-1
Тогда уравнение (11) можно аппроксимировать следующей системой линейных алгебраических уравнений (СЛАУ) вида
Я,
■Ф '
к а --
Я,
(19)
где I - единичная матрица размером ^ _ 1) х (N _ 1). Решая эту систему, получаем вектор ка , проекции которого являются оценками для значений к (¿г-), 1 -1... N _ 1, идентифицируемой ИПФ системы.
Следует сказать, что СЛАУ (19) имеет хорошо обусловленную матрицу (число обусловленности матрицы системы не превышает нескольких десятков), и точность идентификации определяется только ошибками дифференцирования входного (элементы матрицы Ф ') и выходного (проекции вектора / ') сигналов. Вычисление параметра сглаживания в используемых СКС на основе критерия оптимальности (приводимого ниже) дает определенную уверенность в минимальных ошибках вычисления производных от зашум-ленных входного и выходного сигналов идентифицируемой системы, а следовательно, на наименьшую ошибку идентификации ИПФ в целом.
Покажем работу изложенного алгоритма идентификации на примере оценивания ИПФ ПС2 «Вентилятор». Входным сигналом этой подсистемы в процедуре идентификации является выходной сигнал /1н(0 ПС1 (см. рис. 2),
а регистрируются значения /1н - /1я (?,) + ) , ' - 1. N . Однако эти значения были сглажены СКС а^ (^) при идентификации ПС1, и эти сглаженные данные р>2. - (^), ' -1.N, (см. рис. 3 - сплошная кривая)
будут рассматриваться как значения входного сигнала при идентификации ПС2 рассматриваемой сложной системы (см. рис. 1). Первая производная
сплайна Я
Ля ,аш
(О участвует в формировании матрицы Ф ' СЛАУ (18). Ре-
гистрируемые значения /2. - ) + ), ' -1. N выходного сигнала
ПС2 (см. рис. 2) также содержат шумы измерения, и для их фильтрации строится СКС $/ а^ (^). По этому сплайну вычисляется вектор первых производных /2, который входит в правую часть СЛАУ (18). Решая систему (18), определяем оценку к2а для ИПФ ПС2. На рис. 5 показаны значения двух оценок ИПФ: точечной кривой - ИПФ, когда вектор производных /2 вычис-
лялся по интерполяционному сплайну (т. е. без сглаживания выходного сигнала); сплошной кривой - ИПФ, построенная по сглаживающему сплайну. Видно, что вторая ИПФ уже не содержит шумовых осцилляций и в полной мере отвечает качественным априорным представлениям о форме и гладкости ИПФ ПС2 «Вентилятор».
Рис. 5. Оценки ИПФ ПС2
Fig. 5. Impulse response estimates in subsystem 2
4. АЛГОРИТМ ВЫБОРА ПАРАМЕТРА СГЛАЖИВАНИЯ СКС
В работе [11] было показана эффективность использования СКС как способ фильтрации шумов различной статистической природы. Однако выбор параметра сглаживания остается решающим моментом в приведенных алгоритмах идентификации. В работе [13] для сглаживания зашумленных данных были рассмотрены несколько алгоритмов выбора параметра сглаживания. Показано, что эффективной (наилучшей) оценкой для а^ является
величина аф, вычисленная на основе статистического критерия оптимальности фильтрации зашумленных данных. В работе [11] было показано, что величина параметра аф может также использоваться и при вычислении первой производной от экспериментальных данных, искаженных шумами различной статистической природы. Обоснование критерия оптимальности и его свойства подробно изложены в работах [6, 13]. Здесь приведем только конечные соотношения, необходимые для вычисления оценки аф .
Введем в рассмотрение статистику
1 N _
РФ (а) = — ^ , (20)
¿=1
где еI = /1 - Б^,а (/7-) - невязка 7-го измерения. Доказано, что если рш (а) при некотором значении а удовлетворяет неравенству
ир <рш(а) <и в , (21)
-,М 1--N
2
то такое значение можно принять в качестве оценки для аор1 (обозначим это значение как аш). В неравенстве (21) величины и в , и в - квантили
22
X2-распределения с N степенями свободы уровней —, 1 -— соответственно.
2 2
Величина — определяет вероятность ошибки первого рода при проверке гипотезы об оптимальности оценки а^ и, как правило, — = 0,05 . Если N > 30,
то квантили и— , и — для В = 0,05 можно приближенно вычислять по
—N 1—N 22
формулам
= N - \,96^Ш , и0 975^ = N + \96-j2N .
и0,025^ = Я - , и0,975,N
Заметим, что вычисление а^ сводится к решению нелинейного уравнения
Рш (а) = N (21)
итерационными алгоритмами. В качестве а^ принимается очередное приближенное решение а(и), которое удовлетворяет неравенству (21). Даже «медленные» итерационные алгоритмы (например, метод дихотомии) вычисляли значение а^ не более чем за 5 итераций. Выполненные исследования точности оценки а^ [13] показали, что сплайн, построенный при а = а^, имеет:
а) ошибку сглаживания, незначительно (на 5... 8 %) превышающую ошибку сглаживания при параметре а = а0рГ (который можно найти только в
вычислительном эксперименте);
б) ошибку сглаживания, значительно (на 15.35 %) меньше по сравнению с другими способами выбора параметра (подробнее см. [9, 10, 13]).
Всё это позволяет сделать вывод о целесообразности использования сглаживающего сплайна с а = а^ для устойчивого вычисления первых производных входного и выходного сигналов идентифицируемой системы. В качестве подтверждения этого вывода приведем (рис. 6) следующие графики из работы [11]: сплошная кривая - относительная ошибка сглаживания данных, искаженных равномерным шумом с относительным уровнем 5^ = 0,10 ; штриховая кривая - относительная ошибка дифференцирова-
ния, штрих-точечная кривая - значения р^ (а); точечные прямые - квантили и0 025180 , и0 975 180 • Значения а , для которых рw (а) находится между точечными прямыми (т. е. удовлетворяется неравенство (20)), могут быть приняты в качестве оценки а^ для оптимального параметра сглаживания, и эти значения находятся в области минимальной ошибки сглаживания и дифференцирования.
Рис. 6. Свойства оптимальности параметра сглаживания <xw
Fig. 6. The optimality of the smoothing parameter aW
Изложенный алгоритм выбора параметра сглаживания требует досто-
2
верного задания дисперсии шума измерения (см. (19)). В реальном эксперименте такая априорная информация может отсутствовать, и в этом случае можно использовать эффективную оценку дисперсии, вычисляемую на основе дискретного преобразования Фурье зашумленного сигнала и подробно исследованную в работе [15].
ЗАКЛЮЧЕНИЕ
Непараметрическая идентификация сложных динамических систем является весьма сложной для практики задачей. С одной стороны, это обусловлено некорректностью задачи решения интегрального уравнения первого рода, с другой стороны, разными типами входных и выходных сигналов идентифицируемых систем. Предлагаемые в работе алгоритмы идентификации позволяют в полной мере учитывать эти особенности и осуществлять устойчивую непараметрическую идентификацию как отдельных составляющих, так и всей системы в целом. Используемый при этом аппарат сглаживающих кубических сплайнов дает возможность построить эффективный алгоритм фильтрации шумов измерений входного и выходного сигналов системы различной статистической природы. Приводимые результаты решения практи-
ческой задачи идентификации сложной системы «воздухонагреватель-вентилятор-помещение» позволяют рекомендовать изложенные алгоритмы для идентификации импульсных переходных функций других сложных стационарных систем.
СПИСОК ЛИТЕРАТУРЫ
1. Мансуров Р.Ш., Рудяк В.Я. Переходные процессы в системе нагреватель-вентилятор при изменении режима работы вентилятора // Известия высших учебных заведений. Строительство. - 2019. - № 3. - С. 50-63. - DOI: 10.32683/0536-1052-2019-723-3-50-63.
2. Мансуров Р.Ш., РудякВ.Я. Экспериментальное изучение процессов теплообмена при переменных режимах работы системы воздухонагреватель-вентилятор // XXXV Сибирский теплофизический семинар : тезисы докладов Всероссийской конференции с элементами научной школы для молодых ученых. - Новосибирск, 2019. - С. 216.
3. Мансуров Р.Ш., Рудяк В.Я. Экспериментальное изучение переходных процессов в системе нагреватель-вентилятор-помещение // Известия высших учебных заведений. Строительство. - 2018. - № 10. - С. 37-50.
4. Сидоров Д.Н. Методы анализа интегральных динамических моделей: теория и приложения. - Иркутск: Изд-во ИГУ, 2013. - 293 с.
5. ТихоновА.Н., АрсенинВ.Я. Методы решения некорректных задач. - М.: Наука, 1986. -
285 с.
6. Воскобойников Ю.Е. Устойчивые алгоритмы непараметрической идентификации динамических систем. - Новосибирск: НГАСУ (Сибстрин), 2019. - 160 с.
7. Воскобойников Ю.Е. Устойчивые алгоритмы решения обратных измерительных задач. - Новосибирск: НГАСУ (Сибстрин), 2007. - 184 с.
8. ЗавьяловЮ.С., КвасовБ.И., МирошниченкоВ.Л. Методы сплайн-функций. - М.: Наука, 1980. - 345 с.
9. Wang Y. Smoothing splines: methods and applications. - Boca Raton, FL: CRC Press, 2011. - 347 p. - (Monographs on Statistics and Applied Probability; vol. 121).
10. Wahba G. Smoothing noisy data with spline functions system // Numerische Mathematik. - 1975. - Vol. 24, N 5. - P. 383-393.
11. Воскобойников Ю.Е., Боева В.А. Исследования эффективности использования сглаживающих кубических сплайнов в задачах непараметрической идентификации // Автоматика и программная инженерия. - 2019. - № 4 (30). - С. 56-64.
12. БалкП.И., ДолгальА.С. Сплайн-сглаживание экспериментальных данных при нулевом медианном значении помех // Автоматика и телемеханика. - 2017. - № 6. - С. 138-156.
13. ВоскобойниковЮ.Е., Преображенский Н.Г., Седельников А.И. Математическая обработка эксперимента в молекулярной газодинамике. - Новосибирск : Наука, 1984.- 238 с.
14. Воскобойников Ю.Е., Боева В.А. Новый устойчивый алгоритм непараметрической идентификации технических систем // Современные наукоемкие технологии. - 2019. - № 5. -С. 25-29.
15. ВоскобойниковЮ.Е., КрысовД.А. Оценивание характеристик шума измерения в модели «сигнал+шум» // Автоматика и программная инженерия. - 2018. - № 3 (25). - С. 54-61.
Воскобойников Юрий Евгеньевич, выпускник кафедры автоматики НГТУ (НЭТИ), доктор физико-математических наук, заведующий кафедрой прикладной математики Новосибирского государственного архитектурно-строительного университета (Сибстрин), профессор кафедры автоматики Новосибирского государственного технического университета. Заслуженный работник Высшей школы РФ. Основное направление научных исследований - методы и алгоритмы решения некорректных задач интерпретации экспериментальных данных, методы и алгоритмы фильтрации сигналов и изображений.
Является автором более 280 научных работ, в том числе пяти научных монографий по решению некорректно поставленных задач. Е-mail: voscob@mail.ru.
Боева Василиса Андреевна, выпускница кафедры автоматики НГТУ (НЭТИ), аспирант Новосибирского государственного архитектурно-строительного университета (Сиб-стрин) по профилю «Системный анализ, управление и обработка информации», ассистент кафедры прикладной математики Новосибирского государственного архитектурно-строительного университета (Сибстрин). Основное направление научных исследований -методы и алгоритмы непараметрической идентификации реальных технических систем. Является автором 18 научных статей. Е-mail: v.boyeva@sibstrin.ru.
Voskoboynikov Yuri E., graduate of the Department of Automation of NSTU (NETI), Doctor of Physics and Mathematics, Head of the Department of Applied Mathematics in the Novosibirsk State University of Architecture and Civil Engineering (Sibstrin), professor at the Department of Automation in the Novosibirsk State Technical University. Honored Worker of the Higher School of the Russian Federation. The main area of his scientific research is methods and algorithms for solving ill-posed problems of interpreting experimental data; methods and algorithms for filtering signals and images. He is the author of over 280 scientific papers, including 5 scientific monographs on solving ill-posed problems. E-mail: voscob@mail.ru.
Boeva Vasilisa A, post-graduate student in "System analysis, control and data interpretation"; teaching assistant, Department of Applied Mathematics, Novosibirsk State University of Architecture and Civil Engineering (Sibstrin). Her research interests are currently focused on non-parametric identification methods and algorithms for real engineering systems. She has 18 scientific articles. Е-mail: v.boyeva@sibstrin.ru.
DOI: 10.17212/1814-1196-2020-4-47-64 Non-parametric identification algorithms for complex engineering systems
Yu.E. VOSKOBOYNIKO V1,2,a, V.A. BOEVA1,b
1 Novosibirsk State University of Architecture and Civil Engineering (Sibstrin), 113 Leningradskaja Street, Novosibirsk, 630008, Russian Federation
2 Novosibirsk State Technical University, 20 K. Marx Prospekt, Novosibirsk, 630073, Russian Federation
a voscob@mail.ru b v.boyeva@sibstrin.ru
Abstract
In a practice, it often happens that complex engineering systems consist of several interconnected different-type simpler subsystems. An adequate model formulation for every subsystem is impractical due to the complexity of physical processes proceeding in the subsystem. In such cases, a non-detailed black-box model is commonly used. For stationary linear systems (or subsystems), the connection between an input and an output of the black-box is defined by the Volterra integral equation of the first kind with an undetermined difference kernel also known as an impulse response in the automatic control theory. It is necessary to evaluate the unknown impulse response to use the black-box model .This statement is a non-parametric identification problem. For complex systems, the problem needs to be solved both for a whole system and for every isolated subsystem that makes identification substantially complex. Formally, impulse response evaluation is a solution of the integral equation of the first kind for its kernel over registered noise-contaminated discrete input and output values. This problem is ill-
*
Received 14 August 2020. Acknowledgments: The reported study was funded by RFBR, project number 20-38-90041.
posed because of possible solution instability regarding measurement noises in initial data. To find a unique stable solution regularizing algorithms are used, but specific input and output signals in impulse response identification experiments do not allow applying computational methods of these algorithms (system of linear equations or discrete Fourier transformation). In this paper, the authors propose two specific-considering identification algorithms for complex engineering systems. In these algorithms, smoothing cubic splines are used for stable calculation of first derivatives of identified system signals. The results of the complex "Heater-Blower-Room" system identification prove the efficiency of algorithms proposed.
Keywords: non-parametric identification problem, integral Volterra equation of the first kind, integral Volterra equation of the second kind, ill-posed problems, smoothing cubic splines, stable identification algorithms, optimal smoothing parameter evaluation, solution of the practical complex system identification problem
REFERENCES
1. Mansurov R.Sh., Rudyak V.Ya. Perekhodnye protsessy v sisteme nagrevatel'-ventilyator pri izmenenii rezhima raboty ventilyatora [Transient processes in the system the heater-fan when changing the operating mode of the fan]. Izvestiya vysshikh uchebnykh zavedenii. Stroitel'stvo = News of higher educational institutions. Construction, 2019, no. 3, pp. 50-63. DOI: 10.32683/0536-10522019-723-3-50-63.
2. Mansurov R.Sh., Rudyak V.Ya. [Experimental study of heat exchanging processes when system heater-fan changing operation]. ХХХV Sibirskii teplofizicheskii seminar : tezisy dokladov Vse-rossiiskoi konferentsii s elementami nauchnoi shkoly dlya molodykh uchenykh [Proceedings All-Russian conference "XXXV Siberian thermophysical seminar"]. Novosibirsk, 2019, p. 216. (In Russian).
3. Mansurov R.Sh., Rudyak V.Ya. Eksperimental'noe izuchenie perekhodnykh protsessov v sisteme nagrevatel'-ventilyator-pomeshchenie [Experimental study of transient processes in the system the heater-fan-room]. Izvestiya vysshikh uchebnykh zavedenii. Stroitel'stvo = News of higher educational institutions. Construction, 2018, no. 10, pp. 37-50.
4. Sidorov D.N. Metody analiza integral'nykh dinamicheskikh modelei: teoriya i prilozheniya [Analysis methods of integral dynamic models: theory and application]. Irkutsk, ISU Publ., 2013. 293 p.
5. Tikhonov A.N., Arsenin V.Ya. Metody resheniya nekorrektnykh zadach [Methods of ill-posed problems solution]. Moscow, Nauka Publ., 1986. 285 p.
6. Voskoboynikov Yu.E. Ustoichivye algoritmy neparametricheskoi identifikatsii dinamiche-skikh sistem [Stable nonparametric identification algorithms of dynamic systems]. Novosibirsk, NGASU (Sibstrin) Publ., 2019. 160 p.
7. Voskoboynikov Yu.E. Ustoichivye algoritmy resheniya obratnykh izmeritel'nykh zadach [Stable algorithms for inverse problems solution]. Novosibirsk, NGASU (Sibstrin) Publ., 2007. 184 p.
8. Zav'yalov Yu.S., Kvasov B.I., Miroshnichenko V.L. Metody splain-funktsii [Methods of splines]. Moscow, Nauka Publ., 1980. 345 p.
9. Wang Y. Smoothing splines: methods and applications. Boca Raton, FL, CRC Press, 2011.
347 p.
10. Wahba G. Smoothing noisy data with spline functions system. Numerische Mathematik, 1975, vol. 24, no. 5, pp. 383-393.
11. Voskoboynikov Yu.E., Boeva V.A. Issledovaniya effektivnosti ispol'zovaniya sglazhivay-ushchikh kubicheskikh splainov v zadachakh neparametricheskoi identifikatsii [Researches of efficiency of using smoothing cubic splines on nonparametric identification problems]. Avtomatika i programmnaya inzheneriya = Automatics and Software Enginery, 2019, no. 4 (30), pp. 56-64.
12. Balk P.I., Dolgal' A.S. Spline smoothing for experimental data under zero median of the noise. Automation and Remote Control, 2017, no. 6 (78), pp. 1072-1086. DOI: 10.1134/ S000511791706008X. Translated from Avtomatika i telemekhanika, 2017, no. 6, pp. 138-156.
13. Voskoboynikov Yu^., Preobrazhensky N.G., Sedelnikov A.I. Matematicheskaja obrabotka jeksperimenta v molekuljarnoj gazodinamike [Mathematical experiment proceeding in molecular gas dynamics]. Novosibirsk, Science Publ., 1984, 238 p.
14. Voskoboynikov Yu^., Boeva V.A. Novyi ustoichivyi algoritm neparametricheskoi identif-ikatsii tekhnicheskikh sistem [New stable algorithm of nonparametric identification of technical systems]. Sovremennye naukoemkie tekhnologii = Modern high technologies, 2019, no. 5, pp. 25-29.
15. Voskoboinikov Yu£., Krysov D.A. Otsenivanie kharakteristik shuma izmereniya v modeli "signal+shum" [Estimation of the noise measurement characteristics in the model "signal+noise"]. Avtomatika i programmnaya inzheneriya = Automatics and Software Enginery, 2018, no. 3 (25), pp. 54-61.
Для цитирования:
Воскобойников Ю.Е., Боева В.А. Алгоритмы непараметрической идентификации сложных технических систем // Научный вестник НГТУ. - 2020. - № 4 (80). - С. 47-64. - DOI: 10.17212/1814-1196-2020-4-47-64.
For citation:
Voskoboynikov Yu.E., Boeva V.A. Algoritmy neparametricheskoi identifikatsii slozhnykh tekhnicheskikh sistem [Non-parametric identification algorithms for complex engineering systems]. Nauchnyi vestnik Novosibirskogo gosudarstvennogo tekhnicheskogo universiteta = Science bulletin of the Novosibirsk state technical university, 2020, no. 4 (80), pp. 47-64. DOI: 10.17212/1814-11962020-4-47-64.
ISSN1814-1196, http://journals.nstu.ru/vestnik Science Bulletin of the NSTU Vol. 80, No 4, 2020, pp. 47-64