Цифровые субоптимальные алгоритмы обработки данных аэродинамического тензометрического эксперимента Текст научной статьи по специальности «Физика»

Том XII

УЧЕНЫЕ ЗАПИСКИ Ц А Г И

19 8 1

№ 2

УДК 533.6.071.087

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

ЭКСПЕРИМЕНТА

Н. А. Колганов, Я• Ш. Флаксман

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

При проведении аэродинамического тензометрического эксперимента одной из основных проблем является проблема информативности тензометричес-ких измерений. Повышение информативности весового эксперимента может быть достигнуто благодаря применению метода максимального правдоподобия [1, 2], который позволяет определять искомые параметры с заданной точностью за минимальное время, или, что то же самое, находить эти параметры с минимальной погрешностью при заданном времени измерения.

Однако, как отмечается в работе [3], даже в случае однокомпонентных аэродинамических тензовесов (весы с невзаимодействующими компонентами) задача оптимального определения величины аэродинамической нагрузки Я0, является многопараметрической и сводится к решению системы нелинейных интегродифференциальных уравнений. Поэтому целесообразным является использование субоптимальных алгоритмов, которые более просты в аппаратурной реализации и имеют достаточно высокий информационный КПД [2]. Описание одного из таких алгоритмов имеется, например, в работе [4].

В данной работе рассматривается несколько субоптимальных алгоритмов обработки сигнала тензовесов при дискретном поступлении информации. Один из рассмотренных алгоритмов, основанный на методе наименьших квадратов (МНК), применим при отсутствии специальных фильтрующих или] демпфирующих устройств и позволяет по начальному участку кривой, описывающей процесс затухания по времени колебаний сигнала тензовесов, восстанавливать значение Яо при неизвестных динамических характеристиках тензовесов (частота собственных колебаний <о и коэффициент затухания X). Другие алгоритмы предполагают использование линейных фильтров или механических демпферов с последующим осреднением выходного сигнала в цифровой форме. (Аналогичная задача в случае осреднения аналоговым путем рассматривается в работе [3]).

Процесс малых колебаний однокомпонентных тензовесов описывается обыкновенным дифференциальным уравнением второго порядка с постоянными коэффициентами:

= у (0) = ;

(1)

здесь ю0 — частота свободных колебаний весов в отсутствие демпфирования, X — коэффициент затухания, к—коэффициент чувствительности тензовесов, определяемый в ходе статической градуировки, /?0 — прикладываемая в момент времени ¿ = 0 постоянная аэродинамическая нагрузка [при t<^0 значение х(()= = 0, при 0 значение х (£) = 1], — неизвестные параметры, описываю-

щие начальное положение тензовесов, (О—аэродинамический шум с нулевым математическим ожиданием.

Согласно [1, 3] широкополосный аэродинамический шум /?ш (£)/£ может быть заменен эквивалентным белым шумом интенсивности Предполагается отсутствие возмущений, приходящих извне через державку тензовесов, поскольку их влияние устранимо с помощью акселерометров [5]. Тогда весовые измерения будут сопровождаться двумя типами шумов: аэродинамическим шумом, превращающимся на выходе тензовесов в аддитивный стационарный нормальный центрированный шум с известной корреляционной функцией [1, 6], и шумом измерительной аппаратуры (электронные шумы усилителя, помехи в соединительных цепях и т. д.), который можно считать стационарным нормальным некоррелированным центрированным шумом.

1. Рассмотрим алгоритм оценивания величины /?0 ПРИ использовании тензовесов без специальных фильтрующих или демпфирующих устройств. В этом случае, с учетом решения уравнения (1), сигнал, поступающий на регистрирующую систему, запишем в виде

где

У (0 = Ч (Т. О + 6 (0.

■ц (Ъ = » + е'и (Ь сое Ы + с з'т со*), а =

(2)

1 Г X ' с= -тг\-

к со

(/?, - я0) + «1.

СО2 = СОд - X2,

"Гт = (Тг. Та. — . 7б) = а> Ь, с, X, со,

7 — вектор неизвестных параметров, ут—транспонированный вектор, £ (¿) — центрированный шум, сопровождающий весовые измерения.

Сигнал у (¿) подвергается квантованию по времени с помощью аналого-цифрового преобразователя (АЦП) достаточно высокой разрядности. Таким образом, определение аэродинамической нагрузки сводится к нахождению пяти неизвестных параметров Ь, с, X и ш по значениям выборки у г = у (£;), г = 1, 2,..., N. Решение этой пятипараметрической задачи осуществляется методом наименьших квадратов с использованием известных алгоритмов нелинейного

Л Л

оценивания [7, 8]. При этом оценки со и X неизвестных параметров со и X, входящих в соотношения (2) нелинейно, рассчитываются методом градиентного (наискорейшего) спуска [8], в процессе реализации которого на каждой итерации (при фиксированных значениях со и X) определение оценок параметров бис, входящих в (2), сводится к решению системы трех линейных алгебраических уравнений с тремя неизвестными.

Как известно [8], сходимость итерационного процесса градиентного спуска зачастую обусловлена в значительной степени выбором начальных приближений параметров со'0' и В данном случае параметр со'0' рассчитывается следующим образом:

т*

(0) _ 2тс дГ - 1 V Ии) ■

здесь т* — число полных периодов колебаний тензовесов в рассматриваемой временной реализации сигнала, — число отсчетов в /-м периоде [между

У-м и (У+1)-м максимальными значениями амплитуды колебаний].

В то же время оценивается по выборке максимальных значений амплитуды колебаний сигнала тензовесов: = у (¿„ах)> ' = 1, 2,..., т* с помощью метода Хартли—Букера [4, 7]. Здесь — моменты времени, в которые амплитуда колебаний достигает своего максимального значения в ¿-м периоде.

-----оЭоОо^ОСоОоОООоОооОоС

_I_____]___I_I

О 0,1 1,0 1,1 1,1 С

Рис. 1

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

Пригодность используемой модели г\ (7, для аппроксимации показаний исследуемых тензовесов проиллюстрирована на рис. 1, где точки соответствуют экспериментальным данным, сплошная кривая — аппроксимация с помощью математической модели (2), параметры которой рассчитаны в результате применения описанного алгоритма при 7*=0,36 с, а штриховой линией обозначено значение соответствующее истинной величине нагрузки МНК позволяет определить дисперсионные характеристики оценок искомых параметров и построить соответствующие доверительные интервалы [7]. Зависимости относительных

Л Л Л Л Л Л

среднеквадратических отклонений а {&}/&, а {Л}/X и о {ш}/а> от Г приведены на рис. 2. Указанные зависимости дают возможность оценить необходимую продолжительность эксперимента для достижения требуемой точности восстановления искомых параметров.

Отметим, что в условиях градуировочного эксперимента, а также при малом уровне аэродинамического шума измерения сопровождаются, в основном, шумом измерительной аппаратуры, который, как уже отмечалось ранее, можно рассматривать как стационарный нормальный некоррелированный шум. В этом случае алгоритм оценивания, основанный на МНК, совпадает с алгоритмом метода максимального правдоподобия (ММП) и тем самым представляет собой оптимальный алгоритм математической обработки [6]. В реальном аэродинамическом эксперименте (при наличии заметного уровня коррелированного аэродинамического шума) изложенный выше алгоритм, являясь, вообще говоря, субоптимальным, позволяет получить оценки динамических характеристик тензовесов (ю и А), которые могут быть использованы в качестве начальных приближений при построении оптимальных алгоритмов ММП.

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

Обозначим через время установления выходного сигнала тензовесов на уровне &. Поскольку для тензовесов типична малая степень демпфирования

и соответственно высокая добротность = то может оказаться

2 Л

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

0.02

\

\

—

через линейный фильтр, представляющий собой цепочку апериодических звеньев первого порядка с одинаковыми постоянными времени Последнее условие, согласно [3], следует из требования минимизации значения Таким образом, передаточная функция я-звенного линейного фильтра имеет вид [3]:

IУп (Р) =

(1 + Г#Р)а

Оптимальный выбор постоянной времени звена Та также определение времени установления подробно анализируются в работе [3]. Ряд соображений по выбору числа звеньев п содержится в [5]. При Ь > ^ математическое ожидание сигнала г(£) на выходе линейного фильтра близко к постоянному значению 8. В этом случае задача определения аэродинамической нагрузки сводится к нахождению оценки математического ожидания выходного сигнала г (¿).

Сигнал на выходе фильтра при может быть представлен в виде

г (*) = » + £ (г1),

где £ (¿) — шум, сопровождающий измерения.

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

щие соотношения для определения дисперсии а^ коэффициента корреляции г$(т) и времени корреляции нормального случайного процесса £ (¿) [3, 6]:

9 (• Яп [• 1\ о

= * 3 С2 (») ^ (-) = -^г ] с- с05 аш> ^=;

О о г

здесь С2(ш) — квадрат частотной характеристики тензовесов и линейного фильтра,

С2(ш)== [(I-«2 ^ + (3] Оптимальные оценки математического ожидания сигнала г(£) на фоне нормально распределенной коррелированной помехи с (¿). полученные ММП, анализируются в работах [6, 9, 10]. Однако в силу достаточной сложности их практической реализации оказывается целесообразным рассмотрение субоптимальных оценок. При этом необходимо указывать степень их неоптимальности, т. е.

Лп

е2{

0.5

0.1

0.05

| 7/7"» -10

50

\

Ч 100

Ю 50 ¡00

Рис. 3

500 N

насколько точность субоптимальных оценок ниже точности оптимальной оценки искомого параметра (в данном случае

Одной из наиболее просто реализуемых субоптимальных оценок математического ожидания при дискретном поступлении информации является выборочное среднее

N

(3)

г=1

где 21 = г ((¡), ¿1 =

1

Т, г=1, 2,..., Л^, N—число экспериментальных

N—1

отсчетов, ¿¡г + Т — интервал наблюдения.

Оценка (3) является несмещенной и состоятельной, а ее дисперсия равна 16, Ю]

-Л7-1

И2

Ь;=0

—гг- I ГI

¿т

N—1

(4)

Представляет интерес сравнение дисперсии (4) с минимально достижимой дисперсией оценки, полученной ММП и определяемой следующим образом [9]:

°ор» = '

1+Г/2Ч

(5)

На рис. 3 сплошными кривыми показана зависимость а2{&у}/с| для фильтра шестого порядка (параметры фильтра, взятые из работы [3], приведены в таблице, <5= 75, к^ (Л) = 0,05). Здесь же штриховой прямой — значение а^ {/а?. Следует отметить наличие минимума в зависимости а- {8-^ )/а| от N при некотором

N = Л^щщ для заданного времени осреднения Т, что согласуется с теоретическими предсказаниями, содержащимися в монографии [10]. Как видно из приведенных на рис. 3 результатов, величина мало отличается от предельно достижимой дисперсии с~р1.

На рис. 4 показаны зависимости и о^п/а? от безразмерного времени

осреднения Т= ЦТ* (сплошные кривые). Значения Л?тщ и при выполне-

нии условий м0 Т*. > 1 и (?>1 слабо зависят от параметров тензовесов ш0 и X. Аналогичным образом могут быть рассмотрены тензометрические системы с механическим демпфером, которые, как показано в [3], обладают информационным КПД по сравнению с системами, использующими линейные фильтры.

Рис. 4

В случае когда с помощью такого демпфера реализуется предельный апериодический режим колебаний тензовесов при эффективном коэффициенте затухания X* = <"о> статистические характеристики аэродинамического шума, преобразованного задемпфированными тензовесами, могут быть рассчитаны с помощью следующих соотношений [3, 6]:

2 1 л/ 2 Вр = __ Л/0 О)0, = - ,

с 4 ш0

Лг(х) = (1 + »0|т|)

Зависимости с2 {(1Л,}/а? от N, соответствующие различным значениям щТ, для случая задемпфированных тензовесов показаны сплошными кривыми на рис. 5. Здесь же штриховой прямой показана величина с20 вычисленная по формуле (5). Поведение и Мтт в зависимости от безразмерного времени

осреднения Т=ш0Т проиллюстрировано на рис. 4 (штриховые кривые).

Рассмотрим тактику планирования тензометрического эксперимента, целью которого является получение оценки Од, с заданной точностью в. Если дисперсия

2

выходного шума или ее оценка сверху известны априори, то, полагая

10—.Ученые записки" № 2 125

min £2

—=--, из зависимостей, аналогичных приведенным на рис. 4, можно опре-

делить минимальное время осреднения Т и минимальное требуемое число отсчетов A^min- В том случае, когда величина о? заранее неизвестна, можно воспользоваться ее несмещенной и состоятельной оценкой [10]:

j * о2 j^}

SN = yV (1 — k3) X ~ кэ = '

1 = 1 5

где о2 {в-д,) определяется соотношением (4).

п ио Т* мо h ш0 д °opt

1 2 000 15 201 2 000 6 000 1,29

6 3,4 59,3 13,5 28,9 1,09

Система с демпфером 11,75 2 4,5 1,11

В качестве меры точности оценки 8Л, можно взять величину соответствующего доверительного интервала. Тогда тактика планирования эксперимента

0,5

0,1

а,as

е,С11 5 10 50 100 500 N

Рис. 5

должна быть построена таким образом, чтобы вероятность выполнения соотношения — &|<;е оказалась равной заданной величине а (выбранному уровню достоверности), т. е.

Я{|а„-»|<в}=а. (6)

Для того чтобы получить уравнение планирования, воспользуемся методом осреднения по слабокоррелированной выборке [11]. Выберем шаг дискретизации Д из условия (Д) < 1 так, что выборку ги г2, ..., г^ можно считать практически некоррелированной. В этом случае соотношение (6) перепишем следующим образом [12]:

ы " ^-(7)

Е2

где ¡а (Ы—1) —квантиль распределения Стъюдента с N—1 степенями свободы [7].

Число отсчетов Ы, необходимое для оценивания аэродинамической нагрузки с заданной точностью, определяется из уравнения (7) с помощью следующей итерационной_процедуры [12]. Вначале по выборке, содержащей Л^ отсчетов, вычисляются и по формуле (7) определяется Затем, если ЛГг

берется дополнительно ЛЛ^ экспериментальных отсчетов и по выборке, содержащей Л^ отсчетов, оцениваются уточненные значения и После

этого по формуле (7) находится новое приближение N3 для N. Итерационный процесс повторяется до тех пор, пока не будет выполняться условие <1 N или Д'г_[_1 будет незначительно отличаться от Л*';.

Сравним точность метода осреднения по слабо коррелированной выборке с точностью оценивания ММП. Учитывая слабую коррелированность выборки так, что членами г^ (г'Д) при />-2 в формуле (4) можно пренебречь и что обычно

N>1, можно записать дисперсию оценки математического ожидания § следующим образом:

а2

"ММ =-^[1 + 2ге(Д)]. (8)

При тех же предположениях дисперсия оценки ММП в соответствии с формулой (5) принимает вид

а2 ~ ^ С9)

Таким образом, используя (8) и (9), получаем

1М . Д

°opt

,, [1+2ле(Д)]. ^k

Из приведенного в таблице примера видно, что погрешность определения математического ожидания 8- по слабо коррелированной выборке случайного процесса г а) не намного превышает погрешность оптимальной оценки, полученной ММП, что может быть, вообще говоря, скомпенсировано увеличением времени осреднения Т. Увеличение Т приводит к довольно быстрому росту точности оценки Фд,. Кроме того, определение 9 по слабо коррелированной выборке позволяет существенно сократить объем вычислительной работы.

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

ЛИТЕРАТУРА

1. Г р о д з о в с к и й Г. Л. Приложение метода максимума прав-" доподобия к задачам оптимизации обработки данных аэродинамического эксперимента. „Ученые записки ЦАГИ", т. VIII, № 3, 1977.

2. ГродзовскийГ. Л. Об информационном КПД аэродинамического эксперимента. „Ученые записки ЦАГИ", т. IX, № 2, 1978.

3. Гродзовский Г. Л. Информационные характеристики аэродинамических тензометрических систем. „Ученые записки ЦАГИ", т. X, № 4, 1979.

4. 3 л е н к о Ю. А. Применение метода максимального правдоподобия и метода Хартли—Букера к обработке сигнала тензовесов. Сб. докладов II■ Всесоюзной конференции по методам аэрофизических исследований, т. 2, СО АН СССР, ИТиПМ, Новосибирск, 1979.

5. Б о г д а н о в В. В., Ромашкин С. Т. Многокомпонентные весовые измерительные системы для импульсной аэродинамической трубы НАГИ. Труды ЦАГИ, вып. 1599, 1974.

6. Л е в и н Б. Р. Теоретические основы статистической радиотехники. М., „Советское радио", 1975.

7. Химмельблау Д. Анализ процессов статистическими методами. М., „Мир", 1973.

8. П о л а к Э. Численные методы оптимизации. М., „Мир", 1974.

9. Свешников А. А. Прикладные методы теории случайных функций. М., „Наука", 1968.

10. В и л е н к и н С. Я. Статистическая обработка результатов исследования случайных функций. М., „Энергия", 1979.

11. Мирский Г. Я. Аппаратурное определение характеристик случайных процессов. М., „Энергия", 1972.

12. 3 а ж и г а е в Л. С., К и т ь я н A.A., Р о м а н и к о в Ю. И. Методы планирования и обработки результатов физического эксперимента. М., Атомиздат, 1978.

Рукопись поступила ЦП 1979 Переработанный вариант поступил 9jIV 1980