CC BY
15
УДК 519.6:621.317
ЧИСЛЕННЫЙ АНАЛИЗ ОБРАТНОЙ ЗАДАЧИ ТЕОРИИ ИЗМЕРЕНИЙ
Е.В. Харитонова
В работе рассматривается модель измерений, построенная на базе анализа обратной многоточечной краевой задачи для обыкновенного дифференциального уравнения, приводящая к исследованию интегрального уравнения Фредгольма первого рода. Регуляризации процедуры получения решения осуществляется методом невязки.
Одной из проблем теории и практики динамических измерений является проблема оперативного оценивания детектируемого сигнала в условиях, когда входной сигнал трудно поддается прямому измерению, а выходной содержит значительную часть динамической погрешности (например, при измерении импульсных или других, быстро меняющихся во времени, сигналов).
Существует несколько различных подходов к решению этой задачи, Среди них анализ ампли-тудно-фазовых характеристик системы, методы коррекции, методы модельного контроля и др.
В настоящей работе рассматривается модель измерений, построенная на базе анализа обратной многоточечной краевой задачи для обыкновенного дифференциального уравнения, приводящая к исследованию интегрального уравнения Фредгольма первого рода. Для регуляризации процедуры получения решения используется классический метод А.Н. Тихонова.
Существенное упрощение собственно процедуры регуляризации решения и оценивания параметра регуляризации достигается за счет использования конечно-элементных конечномерных аппроксимаций.
1. Введение
Задачи, возникающие в теории динамических измерений, могут быть условно разделены на две группы - задачи восстановления сигнала и задачи анализа динамической погрешности. Первая из упомянутых задач - определение входного сигнала, искаженного средствами измерений, -в общей постановке, представляет собой обратную задачу теории измерений, которая может быть сформулирована как задача решения операторного уравнения
А ■ и(0 = *(/)
относительно функции и{1) при неточно заданных операторе А и правой части х(/). Известно
(например [2],[3]), что в таких задачах обратный оператор А~1 как правило неограничен, что приводит к неустойчивости численных процедур решения указанного уравнения.
В подавляющем большинстве случаев метод невязки [2] Тихонова А.Н. (или метод регуляризации, например, метод квазирешений [3] Иванова В.К. и т.п.) дает возможность получения приближенного регуляризованного решения указанной задачи с одновременной оценкой точности получаемого решения.
В важном для приложений случае, когда функции и(1) и х(0 связаны дифференциальным соотношением Цх(0] = и(0 ? Ц-] - дифференциальный оператор, приведенное выше операторное уравнение представляет собой интегральное уравнение первого рода (Вольтерра или Фредгольма - в зависимости от постановки задачи), методы решения которого хорошо изучены в теории (например [2, 4, 9] и цитированная там литература).
Тем не менее, численная реализация того или иного метода регуляризации по-прежнему является довольно тонкой задачей, и эффективность применяемых алгоритмов напрямую связана с избранным способом дискретизации задачи. Априори, дискретизация задачи может быть осуществлена в следующих направлениях:
• дискретизация интегрального уравнения с последующей регуляризацией полученной системы алгебраических уравнений;
• построение минимизирующего функционала с последующей его дискретизацией проекционными методами (например, методом Ритца);
• построение минимизирующего функционала с последующей дискретизацией уравнения Эйлера, описывающего необходимые условия экстремума построенного функционала.
В настоящей работе используется второй подход, отличающийся от общепринятых методов дискретизации использованием кусочно-элементных базисов Лагранжа, что приводит к значительному снижению объема вычислений и эффективным процедурам регуляризации.
2. Постановка задачи
Входной сигнал и(0 первичного преобразователя (датчика) недоступен прямому наблюдению и регистрации и восстанавливается по наблюдениям за выходными показателями х(0 измерительных приборов.
Наблюдаемый сигнал х(0 является решением краевой задачи
ц*) = х(л)(0+о»-1 (О^""1} (0+- - -+СО*'«+«ь(0*(0 = (1)
С/у(х) = /у, 7 = 1,...л где и^х) - линейные в С"а Ь] функционалы.
Пусть щ (О, I -1,... п - фундаментальная система решений однородного уравнения Цх) = 0.
Если выполнено условие АеХ то задача (1) однозначно разрешима [1] и задача
восстановления С/(0 по экспериментальным данным х{1) может быть поставлена как задача решения интегрального уравнения Фредгольма I рода
ъ
|0(*,г)и(гУг = х(0, (2)
а
где (3(7, г) - функция Грина краевой задачи (1), *(0 - «исправленный» наблюдаемый сигнал,
ь
даваемый соотношением 5(0 = *(0 ■* (0 + гДе 0 = ~ интер-
а
поляционный многочлен, ассоциированный с граничными условиями рассматриваемой краевой задачи, существование и единственность которого обеспечивается упомянутым выше условием
иЛ9г)
3. Функция Грина
Пусть Цх) = и линейное дифференциальное уравнение п-то порядка
Цх) = х{п) + (?)х{п~Х) + ... + аг (Ох' + % ($)х -и, (3)
их{х\ гУ2(0,..., ип(х) - линейно-независимые линейные функционалы (непрерывные) в Сп[а,Ь]. Однородной краевой задачей для уравнения (3) будем называть задачу
¿(х) = 0; = 7 = 1,...я. (4)
Если ~ фундаментальная система решений (3), то решение задачи (4) имеет
вид:
АО+ (5)
где постоянные с^} -1,... п определяются соотношениями:
\и1(ф1) + с2и1(ф2) + ... + спи1(фп) = 0
................................................ (6)
схип{(рх) + с2ип(<р2) +... + спип(рп) = О
Очевидно, что если ранг матрицы и = V Л<рк) равен п, то у системы (6) существует толь-
ко тривиальное решение с, = 0. Краевая задача (4) при этом имеет только нулевое решение. В этом случае неоднородная краевая задача
L(x) = u, £/y(x) = OJ = l,..,/i (7)
однозначно разрешима для любой непрерывной правой части u(t\ и справедлива
формула
ь
x(t) = ¡G(t,T)u(T)dT, (8)
а
где G(t,r)~ функция Грина дифференциального оператора (4), определяемая следующими условиями:
1. G(t,r) - непрерывна и непрерывно-дифференцируема по переменной х вплоть до (п - 2) -го порядка включительно.
2. Vre[a,6] G(/,r) G(t,r) обладает непрерывными производными порядков п-\ и п по переменной х в каждом из интервалов [а, г) и (г, . При этом производная (и-1) порядка име-
ет в точке х = т скачок, равный 1 (в общем случае - равный
, где an(t) - коэффициент при
старшей производной в (3)):
dn~]G(t,T)
dt
п~ 1
dn~lG(t,T)
/=г+
э/
«-1
=i.
(9)
t=T-
3. В каждом из интервалов [а, г) и (г, 6] функция Грина £(/,г) является решением краевой задачи (4).
4. Функция Грина определяется единственным образом, если только выполнено условие однозначной разрешимости задачи (4).
Пусть (р\,(Р2->--<Рп ~~ фундаментальная система решений однородного уравнения 1(х) = 0 . Тогда
V* е[а,г): - ¿>,(гМ(0 ; V* е (г,й]: С(*,г) = £ Д(г)й(г) .
1=1 1=1
Условия непрерывности функции Грина и первых ее (п - 2) производных дают:
¿[а1(г)-Д(г)М(г) = 0
/=1 п
Х[аДг)-Д(г)]9>;(г) = 0
Г=1
1[а,(г)-Д(т)]^и-2>(т) = 0
0=1
Из условия на скачок (и -1) производной в точке / = г, получаем:
¿[а1(г)-Д(г)]^л"1>(т) = -1.
Уравнения (10) и (11) образуют систему п линейных уравнений относительно неизвестных разностей а1 (г)- Д (г), / = 1 ,...п, определитель которой
(10)
(П)
Pi (pi
(pi
92
<Pn <Pn
*0
в силу линейной независимости функций ср; (т). Краевые условия (4) примененные к функции Грина дают:
i/y(G(i,T)) = 0, / = Подробнее, в силу линейности функционалов Uj(x)9 получаем:
(12)
где //Дг) =
a^r), t е[а9т)
Положим у/(г) = а1-(г)-Д.(т),/ = 19...л. Значения/Дг) однозначно определяются упомянутой выше системой (10)—(11). Поэтому система (12) является системой п линейных (относительно, например а,-(г)) уравнений (Д. -а,), определитель которой V ~ Uj{(pi) не равен нулю
в силу принятых допущений. Следовательно, она однозначно разрешима, чем и завершается построение функции Грина краевой задачи (4).
Функция Грина будет задаваться соотношением
п
ы
п
(13)
1=1
Рассмотрим уравнение п -го порядка с краевыми условиями:
Ь(х) = х{п) + апЛ (1)х{п~Х) +... + аг (ОАО + % (О* = и и}-(х) = х(^-)9 у = 1
Здесь е[а,й], а = /0 <...<*„ = £ .
Пусть г выбрана так, что хк<т < Тогда
2]
X Ар,(0
Полагая, как и выше ^ = - Д, запишем условия непрерывности (10) - (11):
(14)
G(t,x)
i=i п
Отсюда
ы
<э(г)
(15)
где а>(т) - определитель Вронского фундаментальной системы решений однородного уравнения
(14), со,(т) - определители, получающиеся из а>(т) заменой i -го столбца столбцом
0
ч-Ъ
Краевые условия (12) в рассматриваемом случае запишутся в виде:
<Р\ <Р[
92 <Р2
(Рп <Р'п
-1 - <Р{ГХ)
t-x
п
X а1(т)<рг(Ъ) = О 5 = ],...,*;
1=1 «
.1=1
Учитывая выражения Д через ^ и , приходим к системе
п
X = 0 5 = 1,...,А:; 1=1
и п
Е = Е Г, (*>,(*,) 5 = Л +
¡=1 г=1
(16)
Полагая и(р1,...,р „) =
<Р\(Р2) (Рг^Рг) ■■■ <Рп(Рг) 9\{Рп) Фг (Рп) ••• 0>ЯО>И)
и.
а, (г)
, получаем решение системы (16) в виде
/ = 1
(17)
где
^2 (О
2>, (<>?>, (^+1)
7=1
2>у (*>,('«)
7=1
••• <Рп(Ь)
■■■ <Рп(*к) ■■■ <Рп(*к+1)
/-й столбец
Регуляризация
Как уже отмечалось выше, задача (2) относится к классу задач, решения которых неустойчивы к малым изменениям исходных данных.
Для нахождения ее решения применим метод регуляризации А.Н. Тихонова [2]. При этом задача решения интегрального уравнения заменяется задачей минимизации функционала
и 1,2 и ц2
Ма = \Аи - х\\ + а й , где О - стабилизирующий оператор, в настоящей работе взятый в виде
и _
р(т)и2(т) + д(т)и,2(т) ¿/г,
а - параметр регуляризации, который определяется из условия минимума невязки и зависит от точности задания правой части уравнения (2).
5. Основные расчетные соотношения
Известно [2-5], что переход к дискретному аналогу нахождения регуляризованных приближенных решений уравнения (2) можно осуществить различными способами. В настоящей работе принят следующий метод: вариационная задача метода невязки дискретизируется конечно-элементными аппроксимациями с последующим анализом уравнения Эйлера дискретной задачи. Решения этой задачи с соответствующим образом подобранным параметром регуляризации и принимается за приближенное решение рассматриваемой задачи.
Подробнее, положим
(=0
где Яг{г) ~ базис Лагранжа конечно-элементных аппроксимаций [6, 7]. В дальнейшем ограничимся использованием подпространств кусочно-линейных функций, где /£,(/) задаются соотношениями
МО
О, t < , u t > t
t-t
lZ*i±L 1
J 4~1 '
, t,<t<t[M.
При этом функционал Ma(u) заменится функцией Ма(мя(0) = Ма(м0,м15.„ми), где и" - значения функции un(t) в узлах
bib п ь[ (п
\а
Ы
^ / « ^
7=1
ЕчЧЧг)
dr.
Обозначив через Д(0 интеграл JG(/, гЦ (г, получим
("0 -"л)
4 •■• 4А
4i4) • ■ ■ 4j
Mr
1=1
+
{w0 ... un
( \
Д p%+q$)dT ... J(MA + ?
{(МЛ + ^ЛО^
- lLMM-гZWrf + F;
IJ=0
0 _
My = J[4^ + + )] W = °>n >
■=0
а а
Необходимое условие минимума функционала Ма запишется в виде
дМа Л
—Г = = з=0,1, 2, ...п.
В итоге получаем систему (п +1) линейных уравнений
1=0
в форме
п Ь,
/=0
V*
/
Jr - 2 j Я0[е(/,т)Я5(г)г/г
а
dt
которая и подлежит решению.
6. Результаты счета
Предлагаемый алгоритм решения задачи восстановления сигнала был опробован на модельном [8] примере датчика второго порядка, описываемого уравнением: x(t) + axx(t) + aQx(t) ~ u(t) , где а0 = 10 ООО, ^ =100 .
Постоянная времени датчика принималась равной Т = 0,01, коэффициент демпфирования £ = 0,5. При моделировании на вход датчика подавался импульсный сигнал в виде полуволны синусоиды u{t) - sin 500?. На выходе датчика дополнительно присутствовало приведенное гармоническое одночас-тотное шумовое воздействие v(¿) = 0,05 sin 5000/.
На рисунке приведены результаты восстановления функции u(t) (п — 100).
Работа выполнена при частичной поддержке гранта РФФИ-УРАЛ 04-01-96073.
Литература
1. Наймарк М.А. Линейные дифференциальные операторы. - М.: Наука, 1969.
2. Тихонов А.Н., Арсенин В .Я. Методы решения некорректных задач. - М.: Наука, 1986.
3. Иванов В.К., Васин В.В., Танана В.П. Линейные некорректные задачи и их приложения. -М.: Наука. 1978.
4. Гончарский A.B., Черепащук A.M., Ягола А .Г. Численные методы решения обратных задач астрофизики. - М.: Наука. 1978.
5. Лоусон Ч., Хенсон Р. Численное решение задач метода наименьших квадратов. - М.: Наука. 1986.
6. Деклу Ж. Метод конечных элементов. - М.: Мир, 1976.
7. Зенкевич О. Метод конечных элементов в технике. - М.: Мир, 1975.
8. Грановский В.А. Динамические измерения. - Л.: Энергоатомиздат, 1984.
9. Верлань А.Ф., Сизиков B.C. Методы решения интегральных уравнений с программами для ЭВМ. - Киев: Наукова думка, 1978.
Поступила в редакцию 13 декабря 2004 г.
