МИНИСТЕРСТВО НАУКИ И ВЫСШЕГО ОБРАЗОВАНИЯ РОССИЙСКОЙ ФЕДЕРАЦИИ
ФЕДЕРАЛЬНОЕ ГОСУДАРСТВЕННОЕ БЮДЖЕТНОЕ ОБРАЗОВАТЕЛЬНОЕ УЧРЕЖДЕНИЕ
ВЫСШЕГО ОБРАЗОВАНИЯ
«НОВОСИБИРСКИЙ ГОСУДАРСТВЕННЫЙ ТЕХНИЧЕСКИЙ УНИВЕРСИТЕТ»
Автоматики
Кафедра _______________________________________________________________________
(полное название кафедры)
Утверждаю
Жмудь В.А.
Зав. кафедрой _______________
_____________________________
(подпись, инициалы, фамилия)
1 г.
«___» _______________ 202__
ВЫПУСКНАЯ КВАЛИФИКАЦИОННАЯ РАБОТА БАКАЛАВРА
Мизюканова Анна Александровна
_______________________________________________________________________________
(фамилия, имя, отчество студента – автора работы)
Применение операции интегрирования по частям при формировании алгебраической
_______________________________________________________________________________
(тема работы)
системы уравнений в параметрической идентификации
_______________________________________________________________________________
Факультет автоматики и вычислительной техники
_______________________________________________________________________________
(полное название факультета)
27.03.04 «Управление в технических системах»
Направление подготовки _________________________________________________________
(код и наименование направления подготовки бакалавра)
_______________________________________________________________________________
Руководитель
от НГТУ
Автор выпускной
квалификационной работы
Чикильдин Г.П
______________________________________
(фамилия, имя, отчество)
Мизюканова А.А.
______________________________________
(фамилия, имя, отчество)
к.т.н., доцент
______________________________________
(ученая степень, ученое звание)
АВТФ, АА-77
______________________________________
(факультет, группа)
______________________________________
(подпись, дата)
______________________________________
(подпись, дата)
Консультанты по разделам:
________________________________________________________________________________
(краткое наименование раздела)
__________________________________________________________
(подпись, дата, инициалы, фамилия)
________________________________________________________________________________
(краткое наименование раздела)
__________________________________________________________
(подпись, дата, инициалы, фамилия)
________________________________________________________________________________
(краткое наименование раздела)
Новосибирск
________________________________________________________________________________
(краткое наименование раздела)
__________________________________________________________
(подпись, дата, инициалы, фамилия)
1
202__
__________________________________________________________
(подпись, дата, инициалы, фамилия)
МИНИСТЕРСТВО НАУКИ И ВЫСШЕГО ОБРАЗОВАНИЯ РОССИЙСКОЙ ФЕДЕРАЦИИ
ФЕДЕРАЛЬНОЕ ГОСУДАРСТВЕННОЕ БЮДЖЕТНОЕ ОБРАЗОВАТЕЛЬНОЕ УЧРЕЖДЕНИЕ
ВЫСШЕГО ОБРАЗОВАНИЯ
«НОВОСИБИРСКИЙ ГОСУДАРСТВЕННЫЙ ТЕХНИЧЕСКИЙ УНИВЕРСИТЕТ»
Автоматики
Кафедра _______________________________________________________________________
(полное название кафедры)
УТВЕРЖДАЮ
Введите
здесь В.А.
свой текст
Зав. кафедрой __________________
Жмудь
(фамилия, имя, отчество)
__________ 26.02.21
(подпись, дата)
ЗАДАНИЕ
НА ВЫПУСКНУЮ КВАЛИФИКАЦИОННУЮ РАБОТУ БАКАЛАВРА
Мизюкановой Анне Александровне
студенту _________________________________________________________________
(фамилия, имя, отчество)
27.03.04 Управление в технических системах
Направление подготовки ____________________________________________________
(код и наименование направления подготовки бакалавра)
__________________________________________________________________________
Факультет автоматики и вычислительной техники
__________________________________________________________________________
(полное название факультета)
Применение операции интегрирования по частям при формировании
Тема _____________________________________________________________________
(полное название темы выпускной квалификационной работы бакалавра)
алгебраической системы уравнений в параметрической идентификации
__________________________________________________________________________
Исследование особенностей алгоритма
Исходные данные (или цель работы) __________________________________________
параметрической идентификации, детальный анализ свойств функций,
__________________________________________________________________________
формирующих алгебраическую систему уравнений, и исследование параметров
__________________________________________________________________________
этих функций для установления рекомендаций по их выбору
__________________________________________________________________________
__________________________________________________________________________
__________________________________________________________________________
Структурные части работы __________________________________________________
Введение
текст
1. Особенности использования операции интегрирования по частям в задаче
__________________________________________________________________________
__________________________________________________________________________
идентификации
__________________________________________________________________________
2. Функции, формирующие систему уравнений
__________________________________________________________________________
3. Разработка алгоритмов и ПО реализации формирующей функции и ее
__________________________________________________________________________
производных
__________________________________________________________________________
4. Исследование влияния корректирующих параметров на свойства формирующей
текст
функции и ее производных
___________________________________________________________________________
Заключение
___________________________________________________________________________
___________________________________________________________________________
___________________________________________________________________________
___________________________________________________________________________
___________________________________________________________________________
___________________________________________________________________________
___________________________________________________________________________
___________________________________________________________________________
___________________________________________________________________________
___________________________________________________________________________
___________________________________________________________________________
Задание согласовано и принято к исполнению.
Руководитель
от НГТУ
Студент
Чикильдин Г.П
______________________________________
(фамилия, имя, отчество)
Мизюканова А.А.
______________________________________
(фамилия, имя, отчество)
к.т.н., доцент
______________________________________
(ученая степень, ученое звание)
АВТФ, АА-77
______________________________________
(факультет, группа)
______________________________________
(подпись, дата)
______________________________________
(подпись, дата)
Тема утверждена приказом по НГТУ № 900/2
от « 03 »
марта
1 г.
202__
1 г.
изменена приказом по НГТУ № _________ от «____» ___________ 202__
7
ВКР сдана в ГЭК № _______, тема сверена с данными приказа
7
___________________________________________________
(подпись секретаря государственной экзаменационной комиссии по защите ВКР, дата)
_________________________________________________
(фамилия, имя, отчество секретаря государственной
экзаменационной комиссии по защите ВКР)
АННОТАЦИЯ
Пояснительная записка содержит 52 страниц текста, 20 рисунков, 6 таблиц, 10 использованных источников.
ПАРАМЕТРИЧЕСКАЯ ИДЕНТИФИКАЦИЯ,
ДИФФЕРЕНЦИРОВАНИЕ СИГНАЛОВ, ФОРМИРУЮЩАЯ ФУНКЦИЯ,
ИНТЕГРИРОВАНИЕ ПО ЧАСТЯМ, АЛГЕБРАИЧЕСКАЯ СИСТЕМА
УРАВНЕНИЙ, МЕТОД НАИМЕНЬШИХ КВАДРАТОВ
В работе рассматривается задача оценивания параметров модели линейного динамического объекта по экспериментально полученным данным. Затрагивается проблема численного дифференцирования зашумленных сигналов и
путь к ее решению с помощью операции интегрирования по частям.
Целью работы является исследование особенностей алгоритма параметрической идентификации, детальный анализ свойств функций, формирующих
алгебраическую систему уравнений, и исследование параметров этих функций
для установления рекомендаций по их выбору.
Для эффективного достижения цели было разработано ПО реализации
функции и ее производных, ПО для определения эффективной длительности
функций, а также ПО, позволяющее вычислить интеграл свертки.
ОГЛАВЛЕНИЕ
ВВЕДЕНИЕ ................................................................................................................ 5
1. ОСОБЕННОСТИ ИСПОЛЬЗОВАНИЯ ОПЕРАЦИИ ИНТЕГРИРОВАНИЯ
ПО ЧАСТЯМ В ЗАДАЧЕ ИДЕНТИФИКАЦИИ ............................................... 7
2. ФУНКЦИИ, ФОРМИРУЮЩИЕ СИСТЕМУ УРАВНЕНИЙ ......................... 11
2.1. Требования, предъявляемые к формирующим функциям ....................... 12
2.2. Выбор вида формирующей функции ......................................................... 12
3. РАЗРАБОТКА АЛГОРИТМОВ И ПО РЕАЛИЗАЦИИ ФОРМИРУЮЩЕЙ
ФУНКЦИИ И ЕЕ ПРОИЗВОДНЫХ ................................................................. 14
3.1. Общие сведения ........................................................................................... 14
3.2. Пакет прикладных программ реализации ФФ и ее производных ........... 15
3.3. Определение реализаций формирующей функции и ее производных ... 15
3.4. Определение эффективной длительности формирующей функции ...... 16
3.5. Программный модуль вычисления интеграла свертки ............................ 25
4. ИССЛЕДОВАНИЕ ВЛИЯНИЯ КОРРЕКТИРУЮЩИХ ПАРАМЕТРОВ НА
СВОЙСТВА ФОРМИРУЮЩЕЙ ФУНКЦИИ И ЕЕ ПРОИЗВОДНЫХ ........ 28
4.1. Влияние корректирующих параметров формирующей функции и ее
производных на временные характеристики ................................................... 30
4.2. Влияние корректирующих параметров формирующей функции и ее
производных на спектральные характеристики .............................................. 37
ЗАКЛЮЧЕНИЕ ....................................................................................................... 42
СПИСОК ИСПОЛЬЗОВАННЫХ ИСТОЧНИКОВ ............................................. 43
ПРИЛОЖЕНИЯ ....................................................................................................... 44
ВВЕДЕНИЕ
Постановка многих задач для современных производственных систем
полностью или частично сводится к задаче идентификации. Под идентификацией подразумевается определение математической модели объекта управления, где априорной информацией являются измеренные (как правило, с помехами) входной и выходной сигналы.
В работе рассматривается задача параметрической идентификации динамического объекта, который описывается дифференциальным уравнением. А
задача в свою очередь сводится к оцениванию коэффициентов этого уравнения. Обычно основным подходом для определения коэффициентов модели является метод наименьших квадратов (МНК) [1]. Частое применения МНК привело к появлению различных модификаций этого метода: прямой, компенсационный, обобщенный, рекуррентный МНК. Наряду с МНК, существуют и другие методы идентификации, находящие свое применение на практике, например, метод инструментальной переменой, алгоритм Качмажа и прочее [2-4].
Поскольку в качестве априорной информации используются входной и
выходной сигналы идентифицируемого объекта и их производные, а измерению доступны только сами сигналы, то возникает проблема многократного
дифференцирования измеренных с помехами сигналов, что, вообще говоря,
представляет собой некорректно поставленную задачу. Теоретически существует возможность получения производных сигнала посредством нерекурсивных дифференцирующих фильтров [5,6]. Однако эти фильтры физически нереализуемы и их искусственное использование предполагает наличие полной конечной реализации сигнала, подлежащего дифференцированию, что не всегда
возможно.
В работе для замены операции дифференцирования сигналов, при формировании алгебраической системы уравнений, используется известная в математике операция интегрирования по частям. Многократное применение этой
операции позволяет заменить процедуры явного вычисления производных сигналов аналитическим дифференцированием некоторых линейно независимых
5
формирующих функций. Но для корректности использования операции интегрирования по частям, формирующие функции и их производные (с точки зрения идентификации) должны удовлетворять определенным свойствам, а
именно, на концах интервала интегрирования они должны быть равны нулю, а
их амплитудные спектры должны быть близки к амплитудно-частотным характеристикам полосовых фильтров с разнесенными полосами пропускания, что
дополнительно позволит сглаживать помехи, искажающие измеренные сигналы входа и выхода идентифицируемого объекта.
Таким образом, целью работы, помимо выбора вида формирующей
функции, является исследование свойств этой функции и ее производных для
корректного использования операции интегрирования по частям при формирования алгебраической системы в параметрической идентификации.
Для замены дифференцирования входного и выходного сигналов с применением операции интегрирования по частям ставится задача проработать
конкретный алгоритм реализации этой процедуры, ориентированный на применение современных средств вычислительной техники, исследовать влияние
корректирующих параметров формирующей функции и ее производных на
временные и частотные характеристики.
6
1. ОСОБЕННОСТИ ИСПОЛЬЗОВАНИЯ ОПЕРАЦИИ
ИНТЕГРИРОВАНИЯ ПО ЧАСТЯМ В ЗАДАЧЕ ИДЕНТИФИКАЦИИ
Пусть идентифицируемый объект описывается математической моделью
в виде дифференциального уравнения
y
(l )
l
(t ) aq y
( q 1)
q 1
r 1
(t ) bz x( z 1) (t ) .
(1.1)
z 1
Ставится задача идентификации (оценивания) параметров aq , q 1, l и
bz , z 1, r 1 .
С учетом обозначений
a j u j (t ) aq y (q 1) (t ), q j, j 1, l ,
a j u j (t ) bz [ x( z 1) (t )], z j l , j l 1, n ,
u0 (t ) y(l ) (t ), n l r 1,
где u j (t ) j 1, n и u0 (t ) назовем координатами объекта, дифференциальное
уравнение может быть представлено в виде
a1 u1(t ) a2 u2 (t ) ... an un (t ) u0 (t ),
или
n
a j u j (t ) u0 (t ) .
(1.2)
j 1
В результате имеем одно уравнение с n неизвестными и приходим к
необходимости формирования системы mu n уравнений.
Умножим каждое слагаемое уравнения (1.2) на линейно независимые
формирующие функции (ФФ) wi (t τ) , i 1, mu и проинтегрируем по в
пределах 0, t , получим
n
t
t
a j wi (t τ) u j (τ)dτ wi (t τ) u0 (τ)dτ,
j 1
0
0
или
7
i 1, mu
(1.3)
n
a j uij (t ) u0i (t ),
j 1
i 1, mu .
(1.4)
Выражение (1.4) представляет собой систему mu n линейных алгебраических уравнений относительно n искомых параметров a j , j 1, n , которую
можно записать следующим образом
U(t )a u0 (t ) ,
(1.5)
где U(t ) uij (t ), i 1, mu , j 1, n ,
a a j , j 1, n ,
u0 (t ) u0i (t ), i 1, mu .
Поскольку чаще всего априорной информацией для решения задачи параметрической идентификации являются измеренные с помехами входной
x (t ) x(t ) x(t ) и выходной y (t ) y(t ) y(t ) сигналы идентифицируемого
объекта, то с учетом помех x(t ) и y (t ) , уравнение (1.3) принимает вид
n
t
n
t
j 1
0
j 1
0
a j wi (t τ) u j (τ)dτ a j wi (t τ) u j (τ)dτ
t
t
0
0
(1.6)
wi (t τ) u0 (τ)dτ wi (t τ) u0 (τ)dτ.
Задав длительность интервала интегрирования T , выражение (1.6) преобразуется к виду
n
T
T
j 1
0
0
a j wi (T τ) u j (τ)dτ wi (T τ) u0 (τ)dτ v(T ),
i 1, mu ,
где
n
T
T
j 1
0
0
v(T ) a j wi (T τ) u j (τ)dτ wi (T τ) u0 (τ)dτ .
Переходя к дискретному времени T K 1 t , (r 1)t , получим
при i 1, mu :
8
n
K
K
a j dr wi( K 1r ) u jr t dr wi( K 1r ) u0r t vK ,
j 1
r 1
(1.7)
r 1
где
n
K
K
j 1
r 1
r 1
vK кв a j dr wi ( K 1r ) u jr t d r wi ( K 1r ) u0r t ,
неизмеримая обобщенная помеха, а d r , r 1, K и кв – коэффициенты и погрешность квадратурной формулы.
Тогда система линейных алгебраических уравнений (1.5) с учетом помех
примет вид
Ua u0 v ,
(1.8)
где U uij , i 1, mu , j 1, n ,
a a j , j 1, n ,
u0 u0i , i 1, mu ,
v vi , i 1, mu .
Поскольку информация о производных входного и выходного сигнала
неизвестна, а известны только сами сигналы x(t ) и y (t ) , то задача многократного численного дифференцирования измеренных с помехами сигналов является некорректно поставленной.
В математическом анализе [7, с. 602] известна процедура интегрирования по частям, в которой для вычисления интеграла приходится разбивать
подынтегральное выражение на два множителя, из которых первый дифференцируется, а второй интегрируется. Другими словами, в подынтегральном произведении двух функций, одна из которых производная, а вторая первообразная относительно первой, можно дифференциал поменять местами, и если вторая функция представлена аналитически, то ее производную можно взять
«вручную».
9
Проиллюстрируем эту операцию в нашем случае на примере второй производной сигнала u j () y (2) (), 0, T . Применим двукратное интегрирование по частям. После первого интегрирования имеем
T
wi (T τ) y
(2)
T
T
(τ)dτ wi (T τ) y (τ) wi(1) (T τ) y (1) (τ)dτ .
(1)
0
0
0
Интегрируя вновь по частям второе слагаемое, получим
T
T
T
(1)
(1)
(1)
wi (T τ) y (τ)dτ wi (T τ) y(τ) wi(2) (T τ) y(τ)dτ .
0
0
0
Выражение многократной операции интегрирования по частям в общем
случае для интегралов вида (1.3) выглядит следующим образом
T
wi (T τ) y
( j)
(τ) dτ
0
j 1
T
T
(1.9)
wi(r ) (T τ) y( j 1r ) () 0 wi( j ) (T τ) y(τ) dτ.
r 0
0
Однако, с точки зрения идентификации, данный подход будет корректен
если
wi (r ) (0) wi ( r ) (T ) 0, i 1, mu , r 0, j 1 ,
(1.10)
что обеспечит равенство нулю первого слагаемого в правой части выражения (1.9).
Из вышесказанного можно сделать вывод, что качество алгоритма идентификации во многом зависит от вида функций wi (t ), i 1, mu , формирующих алгебраическую систему уравнений, от длительности интервала интегрирования T , а также от количества формируемых уравнений mu .
10
2. ФУНКЦИИ, ФОРМИРУЮЩИЕ СИСТЕМУ УРАВНЕНИЙ
Из уравнения (1.8) понятно, что для решения системы линейных алгебраических уравнений необходимо, чтобы матрица U была невырожденная,
т.е. необходимо убедиться, чтобы значение ее определителя было отлично от
нуля. Однако, если количество формируемых уравнений mu больше порядка n (так матрица U может быть прямоугольной), то необходимо применить первую трансформацию Гаусса, чтобы матрица стала квадратной. При
этом система преобразуется к виду
UTUa UTu0 UTv .
(2.1)
Невырожденность матрицы UTU будет обеспечиваться за счет линейной независимости ее строк, а это будет иметь место, если координаты объекта
u j (t ) , j 1, n , и ФФ wi (t ), i 1, mu будут линейно независимые. Линейная
независимость координат предполагает наличие во входном сигнале объекта
не менее n 2 гармонических составляющих.
Решение системы (2.1) отыскивается МНК, для чего минимизируется
функционал
m
1 u 2
J
vk .
mu k 1
(2.2)
В результате приходим к решению усеченной системы уравнений
T
UT
Ua o U u0
(2.3)
относительно оценок искомых параметров ao .
Поскольку обеспечить наличие необходимого количества гармоник во
входном сигнале не представляет особого труда, то основная задача формирования алгебраической системы уравнений представляется в выборе линейно
независимых функций, формирующих эту систему.
11
2.1. Требования, предъявляемые к формирующим функциям
Для корректности описанного способа интегрирования по частям необходимо выполнение условия (1.10), что накладывает на ФФ и ее производные
требование равенства нулю их значений на концах интервала интегрирования,
а для этого, разумеется, необходимо, чтобы ФФ имела затухающий характер.
Стоит отметить, что линейная независимость строк матрицы UTU будет выполняться наиболее благоприятным образом, если амплитудно-частотные характеристики (АЧХ) ФФ будут представлять собой АЧХ полосовых
фильтров с разнесенными полосами пропускания в диапазоне частот
0, W эффективной длительности АЧХ идентифицируемого объекта
W ( j) .
Рисунок 2.1 – АЧХ объекта
2.2. Выбор вида формирующей функции
Как уже упоминалось, свойства формируемой системы уравнений во
много зависят от вида ФФ. В этом плане, определенный интерес вызывает
функция вида [8]
w(t ) α( m 1) m! t m e αt , m l , α 0,
12
(2.4)
и удовлетворяющая условию (1.10), а ее ( j 1) -ая производная записывается в
форме
( j 1)
w
(t ) α
( j 1)
j
b ji vi (t ), j 1, l 1 ,
(2.5)
i 1
где b ji (1)( j 1) ( j 1)!
( j 1)! (i 1)! ,
а vi (t ) α( m 2i ) (m 1 i )! t ( m 1i ) e αt .
Можно заметить, что функция (2.4) зависит от двух параметров и m ,
которые являются ее корректирующими параметрами. Эти параметры влияют
на вид и расположение функции во временной и частотной областях. Для исследования влияния корректирующих параметров на свойства ФФ необходимо
проанализировать поведение ФФ и ее производных как во временной области,
так и в частотной, получив спектральные характеристики ФФ и ее производных посредством преобразования Фурье.
13
3. РАЗРАБОТКА АЛГОРИТМОВ И ПО РЕАЛИЗАЦИИ
ФОРМИРУЮЩЕЙ ФУНКЦИИ И ЕЕ ПРОИЗВОДНЫХ
3.1. Общие сведения
В качестве языка программирования был выбран Fortran, поскольку изначально язык был ориентирован на научные и инженерные вычисления. Язык
жёстко стандартизирован, что позволяет написанные на нем программы с легкостью переносить на различные платформы. На Fortran’e написано большое
количество различных математических библиотек для матричной алгебры и
решения систем линейных уравнений, библиотеки для решения дифференциальных уравнений, быстрых преобразований Фурье и прочее. Для языка
Fortran доступно большое количество библиотек, содержащих подпрограммы
решения классических вычислительных задач (например, LAPACK, BLAS), а
также задач организации распределённых вычислений (например, MPI).
Крупным производителем компилятора Fortran для операционной системы Windows до 1997 года была корпорация Microsoft. Дальнейшей разработкой компилятора стала заниматься фирма Intel. На сегодняшний день компилятор Intel Fortran Compiler является наиболее широко используемым, поскольку он позволяет оптимизировать код под платформы Intel IA-32, x86_64 и
IA-64. Компилятор легко интегрируются с популярными сторонними компиляторами, средами разработки и операционными системами.
В качестве среды разработки программного обеспечения (ПО) была выбрана среда от компании Microsoft – Visual Studio 2017. Данная среда имеет
множество инструментов для упрощения процедур написания, отладки и
сборки ПО.
Текущий раздел включает в себя описание пакета прикладных программ
реализации ФФ и ее производных, а также описание разработанных алгоритмов
определения эффективной длительности ФФ. В разделе приводится программный модуль вычисления интеграла свертки и модуль для сравнения интегралов
в задачах параметрической идентификации. Данное ПО разработано для дальнейшей интеграции в уже существующего ПО, приведенное в [8].
14
3.2. Пакет прикладных программ реализации ФФ и ее производных
Пакет прикладных программ FF предназначен для реализации ФФ и ее
производных во временной и частотной областях, для определения эффективной длительности временных и частотных характеристик посредством разработанных алгоритмов определения, а также для коррекции (при необходимости) корректирующих параметров ФФ и ее производных.
В пакете программ FF производится:
вычисление реализаций импульсной характеристики (ИХ) ФФ wk( j 1) ,
j 1, n 1 , k 1, K w согласно формуле (2.5) в соответствии с шагом дис-
кретизации по времени t и запись в файл;
определение эффективной длительности импульсной характеристики
(ИХ) Tw старшей производной ФФ w(n) (t ) ;
вычисление реализаций АЧХ ФФ Wk , k 1, N1 посредством преобразования Фурье в соответствии с шагом дискретизации по частоте и запись в файл;
определение максимальной частоты эффективной длительности
АЧХ WФФ старшей производной ФФ w(n) (t ) ;
коррекция (при необходимости) максимальной частоты эффективной
длительности АЧХ ФФ WФФ по частоте АЧХ объекта W .
Структурно программа FF состоит из двух блоков: FFPARNM и
FFIXYX. В первом блоке вычисляются ИХ и АЧХ ФФ и определяется эффективная длительность, во втором реализации ИХ и АЧХ ФФ записываются в
файлы.
3.3. Определение реализаций формирующей функции и ее
производных
Подпрограмма
SUBROUTINE FFSTNM (T, AL, M, NM, FF)
вычисляет значение j -ой ФФ в момент времени t согласно формуле (2.5).
15
Входными параметрами подпрограммы являются следующие:
T – момент времени t (k 1)t ;
AL и M – корректирующие параметры ФФ и m ;
NM – номер функции j .
Выходным параметром подпрограммы является следующий:
FF – значение j -ой ФФ в заданный момент времени t согласно формуле (2.5).
При функционировании подпрограммы FFSTNM используется вспомогательная подпрограмма:
подпрограмма FACT (вычисление факториала).
Поскольку вычисления факториала нет в стандарте языка Fortran, то для
корректного функционирования подпрограммы FFSTNM была написана подпрограмма FACT, которая получает значение факториала согласно формуле
ниже.
n! 1 2 3
(n 1) n ,
(3.1)
где n – целое число, а n ! – факториал этого числа.
Подпрограмма
SUBROUTINE FACT (N, F)
вычисляет произведение всех натуральных чисел от 1 до N.
Входным параметром является целое число N.
Выходным параметром является:
F – подсчитанное значение произведения всех натуральных чисел в
цикле по i , где i 1, N .
Программные реализации вычисления значения ИХ ФФ в конкретный
момент времени представлена в приложении E, а вычисления факториала
представлена в приложении F.
3.4. Определение эффективной длительности формирующей
функции
Ввиду того, что старшая производная ФФ w(n) (t ) имеет минимальную
16
скорость затухания, представляется целесообразным эффективную длительность ИХ всех производных определять по n -ой производной.
Для этой цели введем относительный уровень усечения w и эффективная длительность ИХ ФФ Tw будет ограничена точкой K w , после которой ам( n)
( n)
плитудные значения ФФ по модулю будут меньше, чем w wm
, где wm
–
максимальное значение старшей производной ФФ. Иными словами, необходимо, чтобы выполнялось следующе условие
( n)
w( n) (t ) w wm
.
(3.2)
Определение величины эффективной длительности ИХ 4-ой производной ФФ проиллюстрировано на рисунке ниже.
Рисунок 3.1 – Временная характеристика 4-ой производной ФФ w(4) (t )
На основе реализации временных характеристик ФФ и ее производных в
пункте 3.3 и рисунка 3.1 можно разработать алгоритм определения эффективной длительности ИХ ФФ Tw старшей производной w(n) (t ) .
Первым этапом алгоритма определения эффективной длительности ИХ
( n)
ФФ Tw находится максимальное значение ФФ wm
. Для этого сравнению под-
17
лежат соседние значения wk( n ) , k 1,2, ... максимальное из которых сохраняется.
( n)
wm
max wk(n) , wk(n)1 , k 1,2, ... .
(3.3)
Поскольку амплитуда первого экстремума ФФ самая большая по абсолютной величине, то сравнение соседних значений ФФ прекращается, как
только ФФ начинает убывать.
Подпрограмма
SUBROUTINE FFMAX (DT, AL, M, NM, WM)
вычисляет максимальное значение ФФ.
Входными параметрами являются следующие:
DT – шаг дискретизации t [c];
AL и M – корректирующие параметры ФФ и m ;
NM – номер старшей производной ФФ n .
Выходным параметром является следующий:
WM – максимальное значение n -ой производной ФФ.
При функционировании подпрограммы FFMAX используются вспомогательные подпрограммы:
подпрограмма FFSTNM (вычисление значения ФФ в конкретный момент времени согласно (2.5));
подпрограмма FACT (вычисление факториала).
Листинг подпрограммы FFMAX представлен в приложении D.
Для заданного уровня усечения w можно получить произведение мно( n)
жителей w и wm
, которое задает границы определения эффективной длитель-
ности Tw на оси ординат (см. рис. 3.1).
Поскольку аналитическое выражение функции (2.5) содержит многочлен
j -ой степени, то согласно свойству многочленов он имеет не более j корней, а
следовательно, и j пересечений с осью абсцисс. Данное свойство позволит пропустить все экстремумы ФФ, начав поиск эффективной длительности Tw после
18
последнего пересечения ФФ с осью абсцисс.
Таким образом, на следующем шаге алгоритма находится точка последнего пересечения ФФ с осью абсцисс. Выполнение этого условия обеспечивается, когда значение ФФ изменит свой знак на противоположный. Если произведение двух соседних значений ФФ меньше нуля, а произведение может быть
отрицательным только если сомножители имеют разные знаки, то значит ФФ
пересекла ось абсцисс. При совпадении одного из двух соседних значений ФФ
с точкой пересечения его значение равно нулю.
Для нахождения точки последнего пересечения ФФ с осью абсцисс производится подсчет всех пересечений. После чего осуществляется сохранение
точки, в которой ФФ последний раз пересекает ось абсцисс. Первообразная ФФ
имеет пересечение с нулем только в начале координат, в связи с чем необходимо
выполнение дополнительного условия для первообразной, которое начнет определять эффективную длительность только после того, когда ФФ начнет убывать.
Заключительным шагом алгоритма является определение эффективной
длительности ИХ Tw по уровню усечения w . Поскольку после n -ого пересечения, начиная с максимального значения следует затухание ФФ до нуля в бесконечности, главным в третьем шаге алгоритма является сравнение текущего
значения затухающей ФФ с допустимой границей определения эффективной
длительности. Если в некоторой точке текущее значение ФФ входит в допусти( n)
мые границы w wm
, то эффективная длительность ограничивается этой точ-
кой.
В подпрограмме
SUBROUNINE FFTWNMN (DT, AL, M, NM, DE, KW, TW)
производится определение эффективной длительности ИХ ФФ согласно описанному выше алгоритму.
Входными параметрами подпрограммы являются следующие:
DT – шаг дискретизации t [c];
19
AL и M – корректирующие параметры ФФ и m ;
NM – номер старшей производной ФФ n ;
DE – относительный уровень усечения w .
Выходным параметром является:
KW (TW) – значение эффективной длительности старшей производной
ФФ в точках K w , в секундах Tw .
Структурно подпрограмму FFTWNM можно разделить на следующие
этапы согласно вышеизложенному алгоритму определения:
( n)
вычисление максимального значения старшей производной ФФ wm
;
нахождение точки n -го пересечения ИХ с осью абсцисс;
определение эффективной длительности ИХ ФФ Tw K w по уровню
усечения w .
При функционировании подпрограммы FFTWNM используются вспомогательные подпрограммы:
подпрограмма FFMAX (вычисление максимального значения ФФ);
подпрограмма FFSTNM (вычисление значения ФФ в конкретный момент времени согласно (2.5));
подпрограмма FACT (вычисление факториала).
Вычисление значений ФФ и ее производных в частотной области осуществляется по полученным реализациям временных характеристик ФФ и ее
производных w( j 1) (t ), j 1, n 1 в разделе 3.3, используя прямое преобразование Фурье с алгоритмом быстрого преобразования.
Подпрограммы быстрого преобразования Фурье по Кули N1APFK и
N1AKUL приведены в [9].
Подпрограмма
SUBROUTINE FFSPNM (KW, DT, AL, M, NM, NF, OMM)
вычисляет максимальную частоту эффективной длительности АЧХ старшей
производной ФФ WФФ .
20
Входными параметрами подпрограммы являются следующие:
KW – значение эффективной длительности старшей производной ФФ
K w (в точках);
DT – шаг дискретизации t [c];
AL и M – корректирующие параметры ФФ и m ;
NM – номер старшей ФФ j ;
NF – параметр преобразования Фурье N F 2m , m 1,2,3... .
Выходным параметром подпрограммы является следующий:
OMM – значение максимальной частоты эффективной длительности
АЧХ старшей производной ФФ WФФ .
При функционировании подпрограммы FFSPNM используются вспомогательные подпрограммы:
подпрограмма FFSTNM (вычисление значения ФФ в конкретный момент времени согласно (2.5));
подпрограмма FACT (вычисление факториала);
подпрограмма N1APFK (быстрое преобразование Фурье по Кули);
подпрограмма N1AKUL (алгоритм Кули).
Программная реализация подпрограммы FFSPNM находится в приложении H.
Все подпрограммы для определения эффективной длительности вызываются в подпрограмме FFPARNM.
В подпрограмме
SUBROUNINE FFPARNM (NM, AL, M, DT, KW, TW)
производится вычисление реализаций ФФ во временной и частотной областях,
определение эффективной длительности временных и частотных характеристик посредством разработанных алгоритмов, а также корректирование параметров ФФ при необходимости.
По запросу в подпрограмму вводится следующий параметр:
DE – уровень усечения ИХ ФФ w .
21
Входными параметрами подпрограммы являются следующие:
NM – номер старшей производной ФФ n ;
AL и M – корректирующие параметры ФФ и m ;
DT – шаг дискретизации t [c].
Выходным параметром подпрограммы является:
KW (TW) – значение эффективной длительности старшей производной
ФФ K w (в точках), Tw (в секундах).
В подпрограмме FFPARNM используются следующие уже описанные
вспомогательные подпрограммы:
подпрограмма FFSTNM (вычисление значения ФФ в конкретный момент времени согласно (2.5));
подпрограмма FFMAX (вычисление максимального значения ФФ);
подпрограмма FACT (вычисление факториала);
подпрограмма FFTWNM (определение эффективной длительности ИХ
старшей производной ФФ Tw );
подпрограмма N1APFK (быстрое преобразование Фурье по Кули);
подпрограмма N1AKUL (алгоритм Кули);
подпрограмма FFSPNM (определение максимальной частоты эффективной длительности АЧХ старшей производной ФФ WФФ ).
Все реализации ИХ и АЧХ ФФ и ее производных для удобства записываются в файлы с помощью подпрограммы FFIXYX.
Подпрограмма
SUBROUTINE FFIXYX (NK, AL, M, DT, KN, W, NF, N1, AYX)
позволяет вычислить ИХ и АЧХ ФФ и ее производных, а все значения записать в файлы. В подпрограмме FFIXYX формируются файлы данных
FORT1.DAT и FORT2.DAT.
В FORT1.DAT записываются следующие данные
текущее время t (k 1)t ;
22
значения всех функций w( j 1) (t ), j 1, n 1 в текущий момент времени t .
В FORT2.DAT записываются следующие данные:
текущее значение частоты (k 1) ;
значения спектров всех функций W j ( j) , j 1, n 1 в соответствии с
текущим значением частоты .
Входными параметрами подпрограммы являются следующие:
NK – количество функций w( j 1) (t ), j 1, n 1 ;
AL и M – корректирующие параметры ФФ и m ;
DT – шаг дискретизации t [c];
KN – количество точек, определяемых ИХ ФФ, равное эффективной длительности ФФ K w ;
NF – параметр преобразования Фурье N F 2m , m 1,2,3... ;
N1 – параметр N1 N F 2 1 , определяющий количество дискретных
значений частотных характеристик в диапазоне 0, 0 , где 0 t .
Выходными параметрами подпрограммы являются следующие:
W(KN, NK) – двумерный вещественный массив дискретных значений ИХ
функций w( j 1) (t ), j 1, n 1 ;
AYX(NK, N1) – двумерный вещественный массив дискретных значений
АЧХ функций W j ( j) , j 1, n 1 .
При функционировании программы FFIXYX используются следующие
вспомогательные подпрограммы:
подпрограмма FFSTNM (вычисление значения ФФ в конкретный момент времени согласно (2.5));
подпрограмма FACT (вычисление факториала);
подпрограмма N1APFK (быстрое преобразование Фурье по Кули);
подпрограмма N1AKUL (алгоритм Кули).
23
Весь пакет программ приведен в приложении. Листинг программы FF
приведен в приложении A. Программно-реализованные алгоритмы определения эффективной длительности временных и частотных характеристик старшей производной ФФ вида (2.5) FFPARNM представлены в приложении B. А
подпрограмма вычисления всех значений ФФ и ее производных во временной
и частотной областях, а также запись этих значений в файлы FFIXYX представлена в приложении H.
Графики спектральных характеристик в заданном диапазоне частот
представлены ниже. Временным характеристикам ФФ и ее производных
w( j 1) (t ) соответствуют спектральные характеристики W j ( j) , где
j 1, n 1 .
Рисунок 3.2 – Спектральные характеристики ФФ и ее производных
при 3 , m 5
В зависимости от корректирующих параметров функции вида (2.4) меняются свойства ФФ и ее производных, меняется и расположение спектров в заданном диапазоне частот. Следовательно, есть необходимость исследовать
влияние корректирующих параметров на свойства ФФ во временной и частотной областях.
24
3.5. Программный модуль вычисления интеграла свертки
Поскольку для описанного подхода замены дифференцирования сигналов с помощью операции интегрирования по частям (раздел 1) необходимо вычисление интегралов, то необходимо проработать алгоритм определения интеграла свертки. В подразделе описывается программная реализация интеграла
свертки на конечном интервале времени.
Любой определенный интеграл представляет собой площадь, ограниченную подынтегральной функцией. Вычислить эту площадь можно как сумму
площадей некоторых элементарных геометрических фигур, например, прямоугольников.
Рисунок 3.3 – Вычисление определенных интегралов
Более точный результат интегрирования будет иметь место, если в качестве элементарных геометрических фигур использовать не прямоугольники, а
трапеции.
Интеграл свертки имеет вид
T
y (t ) x() w(T )d .
(3.4)
0
Согласно рисунку 3.3 значение интеграла свертки (3.4) будет определятся как
K x Kw
Y t xi w j ,
i 1 j 1
25
(3.5)
где K x – количество точек функции x(t ) , K w – количество точек функции
w(t ) , а t - шаг дискретизации.
В программе INTEGRAL производится:
вычисление реализаций ИХ ФФ wk( j 1) , j 1, n 1 , k 1, K w согласно формуле (2.5) в соответствии с шагом дискретизации по времени t и
запись в файл;
определение эффективной длительности ИХ старшей производной
ФФ Tw ;
вычисление реализаций ИХ некоторого синусоидального сигнала xk( r ) ,
r 0, n k 1, K x в соответствии с шагом дискретизации по времени t и
запись в файл;
вычисление интегралов свертки F1 и F2 .
По запросу в подпрограмму вводятся следующие параметры:
ALF и MF – корректирующие параметры ФФ и m ;
DE – уровень усечения ИХ ФФ w .
Промежуточными параметрами подпрограммы являются следующие:
NK – количество функций (первообразная и производные ФФ);
DT – шаг дискретизации t [c].
Выходными параметрами подпрограммы являются:
F1(NK) – одномерный массив с вычисленными значениями интеграла
свертки всех производных сигнала x(r ) (t ), r 1, n и первообразной ФФ
w(t ) ;
F1(NK) – одномерный массив с вычисленными значениями интеграла
свертки всех производных ФФ w( j 1) (t ), j 2, n 1 и первообразной сигнала x(t ) .
В программе INTEGRAL используются следующие вспомогательные
подпрограммы:
26
подпрограмма FFSTNM (вычисление значения ФФ в конкретный момент времени согласно (2.5));
подпрограмма FACT (вычисление факториала);
подпрограмма FFMAX (вычисление максимального значения ФФ);
подпрограмма FFTWNM (определение эффективной длительности ИХ
старшей производной ФФ Tw );
подпрограмма N1YSWI (вычисление интеграла свертки на конечном
времени T ).
Подпрограммы FFSTNM и FACT были описаны подразделе 3.3. Подпрограммы FFMAX и FFTWNM в подразделе 3.4.
Подпрограмма
SUBROUTINE N1YSWI (KW, W, KON, X, Y)
вычисляет значение интеграла свертки на конечном интервале времени.
Входными параметрами подпрограммы являются следующие:
KW – значение эффективной длительности K w первой подынтегральной
функции w(t ) (в точках);
W(KW) – массив значений ИХ первой подынтегральной функции w(t ) ;
KON – количество точек второй подынтегральной функции x(t ) ;
X(KON) – массив значений временной характеристики второй подынтегральной функция x(t ) .
Выходным параметром подпрограммы является:
Y – сумма значений произведений подынтегральных функций x(t )
и w(t ) .
В результате значение интеграла вида (3.4) вычисляется, как произведение шага дискретизации t на сумму значений произведений подынтегральных функций x(t ) и w(t ) .
27
4. ИССЛЕДОВАНИЕ ВЛИЯНИЯ КОРРЕКТИРУЮЩИХ
ПАРАМЕТРОВ НА СВОЙСТВА ФОРМИРУЮЩЕЙ ФУНКЦИИ
И ЕЕ ПРОИЗВОДНЫХ
Прежде, чем использовать функцию (2.5) для формирования алгебраической системы в параметрической идентификации, рассмотрим пример корректности описанного подхода интегрирования по частям с выполнением
условий (1.10).
Пусть мы имеем некоторый, измеренный без помех, сигнал x(t ) , аналитическое выражение которого x(t ) sin(2t ) .
Проинтегрируем по частям первую производную входного сигнала
x(1) (t ) 2cos(2t ) и ФФ вида (2.5) на интервале 0, T , получим
(m1)
m (T )
2cos(2)d .
w(T ) x ()d m! (T ) e
T
T
(1)
0
0
Теперь проинтегрируем первую производную ФФ и сигнал x(t ) на том
же интервале интегрирования, получим
T
w
(1)
(T ) x()d
0
( m1)
m
α e(T )
(T )m
(T )( m1) sin(2)d .
(m 1)!
m!
0
T
Согласно выражению (1.9) значения интегралов должны совпасть, если
выполняется условие (1.10). В общем виде эти интегралы можно обозначить
как
T
F1r w(T ) x
(r )
T
()d ,
F2r w(r ) (T ) x()d ,
0
0
где r – номер производной.
Пусть длительность интервала интегрирования T Tx Tw K w 1 t ,
где Tx – длительность сигнала x(t ) , Tw – эффективная длительность ФФ, t –
28
шаг дискретизации по времени.
Выполнение условия (1.10) требует предварительного оценивания эффективной длительности ФФ согласно рисунку 3.1.
Положим r 4 . Зададим t 0.01 с , а корректирующие параметры ФФ
w(r ) (t ), r 0, 4 положим 3 , m 5 .
В таблице 4.1 приведены значения 4-х производных ФФ на концах интервала интегрирования, определяемые соответствующими уровнями усечения w (0.01, 0.001, 0.0001 и 0.00001). Поскольку wi (r ) (0) 0, i [1, mu ] ,
r 0, n априори, то в таблице 4.1 приведены значения лишь для w(r ) (Tw ),
r 0, 4 .
Таблица 4.1 – Влияние w на значения ФФ и ее производных в момент
времени Tw
w
Tw
w(Tw )
w(1) (Tw )
w(2) (Tw )
w(3) (Tw )
w(4) (Tw )
0.01
3.670
0.6688E-01
0.1095E+00
0.1545E+00
0.1582E+00
0.9458E-02
0.001
5.590
0.1728E-02
0.3638E-02
0.7383E-02
0.1428E-01
0.2585E-01
0.0001
6.840
0.1115E-03
0.2529E-03
0.5619E-03
0.1217E-02
0.2557E-02
0.00001
7.930
0.8873E-05
0.2102E-04
0.4911E-04
0.1128E-03
0.2543E-03
Таблица 4.1 демонстрирует, как уменьшение уровня усечения w , приводит к росту эффективной длительности ИХ ФФ Tw , тем самым улучшая выполнение условия (1.10).
Теперь согласно эффективной длительности 4-ой производной ФФ Tw из
таблицы 4.1 сравним значения интегралов F1r и F2r r 1, 4 в зависимости
от уровня усечения w для первых четырех производных ФФ w(t ) и сигнала x(t ) .
29
Таблица 4.2 – Результаты сравнения F1r и F2r в зависимости от w
1w 0.01
2w 0.001
3w 0.0001
4w 0.00001
r
F1r
F2r
F1r
F2r
F1r
F2r
F1r
F2r
1
-0.55396
-0.55462
0.13229
0.13227
-0.49581
-0.49581
0.64548
0.64548
2
0.74786
0.88161
-1.30188
-1.29843
0.88225
0.88247
0.30824
0.30826
3
2.21583
1.99946
-0.52914
-0.53635
1.98324
1.98274
-2.58193
-2.58197
4
-2.99143
-3.22268
5.20753
5.21445
-3.52900
-3.53421
-1.23298
-1.23203
Следуя полученным в таблице данным, видно, что чем меньше порядок
производной, тем лучше совпадают значения интегралов. Для первой производной при 3w и 4w интегралы полностью совпадают с точностью до стотысячных, при 2w до десятитысячных, а при 1w до тысячных. Начиная со второй производной, уровня 1w 0.01 оказывается не достаточно для выполнения условия (1.10), поскольку разница заметна уже в первом знаке после запятой. Для третьей производной уровни w 0.01 обеспечивают совпадение интегралов до двух-четырех значащих цифр после запятой.
Вышесказанное дает сделать вывод, что операция замены дифференцирования сигнала аналитически взятой производной ФФ путем интегрирования
по частям выполняется с достаточной точностью (совпадение интегралов до
двух-трех значащих цифры после запятой) при выборе уровня усечения w 0.001 .
4.1. Влияние корректирующих параметров формирующей функции
и ее производных на временные характеристики
Рассмотрим ФФ и ее производные во временной области для дальнейшего исследования влияния корректирующих параметров на свойства ФФ. Рисунки 4.2–4.4 иллюстрируют функцию (2.5) и ее 4 производных на интервале
определения при значениях корректирующих параметров m const 5 ,
и var .
30
Рисунок 4.1 –ФФ и ее производные ( 3 , m 5 )
Рисунок 4.2 – ФФ и ее производные ( 4 , m 5 )
31
Рисунок 4.3 – ФФ и ее производные ( 5 , m 5 )
Рисунок 4.4 – ФФ и ее производные ( 6 , m 5 )
Из графиков отчетливо видно, что чем меньше порядок производной,
тем меньше максимальное значение ее амплитуды и тем выше скорость ее затухания.
Рисунки 4.5–4.8 иллюстрируют временные характеристики ФФ и 4-х ее
производных при значениях корректирующих параметров const 3 ,
и m var .
32
Рисунок 4.5 – ФФ и ее производные ( 3 , m 4 )
Тут стоит отметить, что при α 3 , m 4 четвертая производная ФФ выходит не из начала координат, что нарушает условие (1.10). Поэтому такая комбинация корректирующих параметров ФФ может быть применима только для
объектов не старше 3-го порядка.
Рисунок 4.6 –ФФ и ее производные ( 3 , m 5 )
33
Рисунок 4.7 – ФФ и ее производные ( 3 , m 6 )
Рисунок 4.8 – ФФ и ее производные ( 3 , m 7 )
Предложенные характеристики наводят на мысль о том, что есть необходимость детально исследовать влияние корректирующих параметров на свойства ФФ и ее производных. Рисунки показывают, что увеличение приводит
34
к сжатию временных характеристик вдоль оси абсцисс, что и следовало ожидать, т.к. величина этого параметра определяет скорость затухания ( e αt ).
Увеличение m , наоборот, растягивает временную характеристику функции
вдоль оси абсцисс.
Для выполнения условия (1.10) необходимо оценить влияние корректирующих параметров ФФ ( и m ) на значение функций в конечный момент
времени.
Исследуем влияние на значение ФФ и ее производных в конечный момент времени при уровне усечения w 0.001, шаге дискретизации t 0.01
и m 5 – таблица 4.3, а также сравним значения интегралов F1r и F2r
при r 1, 4 – таблица 4.4.
Таблица 4.3 – Влияние на Tw и w(r ) (Tw ) , r 0, 4
Tw
w(Tw )
w(1) (Tw )
w(2) (Tw )
w(3) (Tw )
w(4) (Tw )
3
5.590
0.1728E-02
0.3638E-02
0.7383E-02
0.1428E-01
0.2585E-01
4
4.200
0.2256E-02
0.6337E-02
0.1717E-01
0.4433E-01
0.1072E+00
5
3.360
0.2820E-02
0.9902E-02
0.3353E-01
0.1082E+00
0.3270E+00
6
2.800
0.3384E-02
0.1426E-01
0.5793E-01
0.2244E+00
0.8138E+00
Приведенная таблица показывает, что использование параметра
для производных 3-го, 4-го порядка требует меньшего уровня усечения и уровень w 0.001 является недостаточным для выполнения условия (1.10).
Таблица 4.4 – Результаты сравнения интегралов F1r и F2r в зависимости от
4
5
6
3
r
F1r
F2r
F1r
F2r
F1r
F2r
F1r
F2r
1
0.13229
0.13227
0.80471
0.80469
-0.34964
-0.34967
-1.26078
-1.26081
2
-1.30188
-1.29843
1.26239
1.26690
2.46473
2.47036
1.46768
1.47445
3
-0.52914
-0.53635
-3.21884
-3.23143
1.39858
1.37889
5.04312
5.01475
4
5.20753
5.21445
-5.04955
-5.06241
-9.85891
-9.86943
-5.87073
-5.53668
35
Из таблиц 4.3-4.4 видно, что для уровня w 0.001 условие (1.10) для
всех 4-х производных ФФ достигается наилучшим образом при наименьшем
значении 3 , а совпадение интегралов F1r и F2r выполняется до сотых.
Схожие результаты при 4 и 5 . При значении параметра больше пяти
ухудшается точность для третьей и четвертой производных и, судя по результатам, в данном примере можно ограничиться значением 3, .
При тех же условиях и оценим влияние параметра m на выполнение соотношения (1.10) – результаты в таблице 4.5.
Таблица 4.5 – Влияние m на Tw и w(r ) (Tw ) , r 0, 4
m
Tw
w(Tw )
w(1) (Tw )
w(2) (Tw )
w(3) (Tw )
w(4) (Tw )
4
3.150
0.7844E-01
0.1357E+00
0.2032E+00
0.2221E+00
-
5
5.590
0.1728E-02
0.3638E-02
0.7383E-02
0.1428E-01
0.2585E-01
6
6.670
0.5458E-03
0.1146E-02
0.2334E-02
0.4572E-02
0.8510E-02
7
7.560
0.2597E-03
0.5385E-03
0.1085E-02
0.2110E-02
0.3923E-02
Таблица 4.6 – Результаты сравнения интегралов F1r и F2r в зависимости от m
m4
m5
m6
m7
r
F1r
F2r
F1r
F2r
F1r
F2r
F1r
F2r
1
-0.81862
-0.81940
0.13229
0.13227
-0.54135
-0.54136
-0.25176
-0.25176
2
0.26007
0.41695
-1.30188
-1.29843
-0.22077
-0.21968
0.76856
0.76908
3
3.27447
3.00613
-0.52914
-0.53635
2.16540
2.16313
1.00702
1.00596
4
-
-
5.20753
5.21445
0.88309
0.88339
-3.07425
-3.07416
В данном случае, чем больше параметр m , тем лучше выполняется соотношение (1.10). При m 7 наблюдается практически полное совпадение значений интегралов F1r и F2r , а результаты по длительности функций иллюстрируют возможность их использования при m 4 .
36
4.2. Влияние корректирующих параметров формирующей функции
и ее производных на спектральные характеристики
Для того, чтобы АЧХ объекта и АЧХ ФФ располагались в одном диапазоне частот, рассмотрим их расположение в зависимости от корректирующих
параметров ФФ, например, в диапазоне частот 0, 20 рад/с.
Рисунки 4.9–4.12 иллюстрируют расположение спектров ФФ и 4-х ее
производных при значениях корректирующих параметров m const 5 ,
и var .
Рисунок 4.9 – Спектральные характеристики ФФ ( 3 , m 5 ) и ее
производных
37
Рисунок 4.10 – Спектральные характеристики ФФ ( 4 , m 5 ) и ее
производных
Рисунок 4.11 – Спектральные характеристики ФФ ( 5 , m 5 ) и ее
производных
38
Рисунок 4.12 – Спектральные характеристики ФФ ( 6 , m 5 ) и ее
производных
По частотным характеристикам, подобно временным, видно, что чем
меньше порядок производной ФФ n , тем быстрее приближается к нулю ее значение.
Рисунки 4.13–4.16 иллюстрируют спектральные характеристики ФФ и 4х ее производных при значениях корректирующих параметров const 3 ,
и m var .
39
Рисунок 4.13 – Спектральные характеристики ФФ ( 3 , m 4 ) и ее
производных
Рисунок 4.14 – Спектральные характеристики ФФ ( 3 , m 5 ) и ее
производных
40
Рисунок 4.15 – Спектральные характеристики ФФ ( 3 , m 6 ) и ее
производных
Рисунок 4.16 – Спектральные характеристики ФФ ( 3 , m 7 ) и ее
производных
По приведенным частотным характеристикам видно, что увеличение параметров и m приводит к уменьшению диапазона расположения спектров
ФФ.
41
ЗАКЛЮЧЕНИЕ
Основной проблемой в решении задачи параметрической идентификации
динамических объектов является измерение производных входного и выходного сигналов идентифицируемого объекта, что представляет собой априорную
информацию для решения этой задачи.
В работе предложен подход, базирующийся на замене дифференцирования зашумленных сигналов посредством использования известной операции
интегрирования по частям при формировании системы линейных алгебраических уравнений относительно идентифицируемых параметров.
С точки зрения параметрической идентификации сформированы требования, на основании которых предложен вид импульсных характеристик фильтров, формирующих систему уравнений.
Разработано ПО, позволяющее определять импульсные характеристики
функций и их производных, используемые при интегрировании по частям.
Разработано ПО, проведены исследования и даны рекомендации по выбору параметров формирующих фильтров, определяющих временные и частотные свойства импульсных характеристик и их производных.
Разработаны алгоритмы и ПО, позволяющие корректировать параметры
формирующих фильтров в процессе формирования алгебраических уравнений
в зависимости от частотных свойств идентифицируемого объекта.
Приведены результаты функционирования, разработанного ПО, в том
числе и в правомерности использования операции интегрирования по частям
при замене численного дифференцирования сигналов аналитическим дифференцированием импульсных характеристик формирующих функций.
42
СПИСОК ИСПОЛЬЗОВАННЫХ ИСТОЧНИКОВ
1. Турчак Л.И. Основы численных методов. - М.: Наука, 1987. – 320 с.
2. Огарков М.А. Методы статистического оценивания параметров случайных процессов. М.: Энергоатомиздат, 1990. – 208 с.
3. Эйкхофф П. Основы идентификации систем управления. – М.: Мир,
1975. – 686 с.
4. Чикильдин Г.П. Идентификация динамических объектов: учеб. пособие. – Новосибирск: Изд-во НГТУ, 2017. – 88 с.
5. Анисимов А.С. Методы цифровой фильтрации: учеб. пособие. – Новосибирск: НЭТИ, 1992 – 82 с.
6. Кононов В. Т., Худяков Д. С., Чикильдин Г. П. Цифровая фильтрация
сигналов. – Новосибирск: Изд-во НГТУ, 2008. – 64 с.
7. Фихтенгольц, Г. М. Курс дифференциального и интегрального исчислен: Том 2 / Г. М. Фихтенгольц. – М.: Наука, 1966. – 800 с.
8. Анисимов А. С. Исследование алгоритмов идентификации импульсной и частотных характеристик – Новосибирск: Изд-во НГТУ, 1996. – 48 с.
9. Анисимов А. С., Кононов В. Т, Чикильдин Г. П. Алгоритмы преобразования линейных динамических моделей: учеб. пособие. – Новосибирск: Издво НГТУ, 2004. – 146 с.
10. Анисимов А. С., Чикильдин Г. П., Худяков Д. С. Методы обработки
сигналов: учеб. пособие. – Новосибирск: Изд-во НГТУ, 2018. – 107 с.
43
ПРИЛОЖЕНИЯ
Приложение A
PROGRAM FF
USE IFNLS
Integer CP
REAL
*
*
*
A(7),B(3),AO(5),W(2049),H(1025),WS(1000),
DS(10000),SX(10000),UX(10000),UY(10000),
WF(2000,5),AYX(5,1025),XS(1000),YS(1000),
U(500,5),UTR(5,500),U0(500),V(7,7),V0(7),E(7)
CP = NLSSetLocale("Russian", "Russia", 1251)
C
Запрос необходимых параметров формирующей функции
PRINT *, 'Порядок ФФ N ->'
READ *, N
NK=N+1
PRINT *, ' ALF ->'
READ *, ALF
PRINT *, ' MF ->'
READ *, MF
DT=0.01
C
C
C
Определение эффективной длительности ИХ и АЧХ
старшей производной формирующей функции
и коррекция ее параметров
CALL FFPARNM (NK,ALF,MF,DT,KWF,TWF)
C
C
Вычисление ИХ и АЧХ формирующей функции и
ее производных и запись в файл
NF=1024
N1=NF/2+1.01
CALL FFIXYX (NK,ALF,MF,DT,KWF,WF,NF,N1,AYX)
STOP
END
Приложение B
SUBROUTINE FFPARNM (NM,AL,MF,DT,KW,TW)
44
PRINT*, 'Уровень усечения ИХ DEL->'
READ* , DE
1
CONTINUE
C
Эффективная длительность TW ИХ ФФ_NM
CALL FFTWNM (DT,AL,MF,NM,DE,KW,TW)
PRINT 101,TW,KW
C
Эффективная длительность OMM SP ФФ_NM
PRINT*,' Параметр БПФ NF->'
READ *, NF
CALL FFSPNM (KW,DT,AL,MF,NM,NF,OMM)
PRINT*,' Макс. частота эффект. длит. АЧХ ФФ_NM
OMM=', OMM
C
Коррекция параметров ФФ
PRINT*,' Изменить параметр MF ФФ_NM ? (1 / 0)'
PRINT*,' Условие изменения параметра MF->'
READ *, IP
IF (IP.EQ.0) GO TO 5
PRINT*, 'Параметр ФФ MF->'
READ* , MF
GO TO 1
5
PRINT*,' Изменить параметр AL ФФ_NM ? (1 / 0)'
PRINT*,' Условие изменения параметра AL->'
READ *, IP
IF (IP.EQ.0) GO TO 6
PRINT*, 'Параметр ФФ AL->'
READ* , AL
GO TO 1
CONTINUE
6
101
FORMAT (2X,'Эффект. длит. ИХ ФФ_NM TW=',F7.3,'
KW=',I4)
RETURN
END
Приложение C
SUBROUTINE FFTWNM (DT,AL,M,NM,DE,KW,TW)
45
C
C
Определение длительности Tw (Kw) ИХ
NM-ой производной формирующей функции на уровне DE
C
Вычисление максимального значения ИХ
CALL FFMAX (DT,AL,M,NM,WM)
DE=WM*DE
C
C
Нахождение точки NM-го пересечения ИХ ФФ
с осью абсцисс
ISY=1
K=0
W1=0.
K=K+1
T=(K-1)*DT
CALL FFSTNM (T,AL,M,NM,FF)
W2=FF
IF (NM.EQ.1.AND.W2.LT.W1) GO TO 4
IF ((W1*W2).GE.0.) GO TO 2
ISY=ISY+1
IF (ISY.GE.NM) GO TO 3
CONTINUE
W1=W2
GO TO 1
1
2
3
KW=K+1
TW=(KW-1)*DT
FLAG=0
C
4
Определение длительности ИХ на уровне DE
K=K+1
W1=W2
T=(K-1)*DT
CALL FFSTNM (T,AL,M,NM,FF)
W2=FF
IF (FLAG.EQ.1) THEN
IF (ABS(W1).LE.DE) THEN
KW=K
TW=(KW-1)*DT
GOTO 5
ELSE
GOTO 4
END IF
END IF
46
IF(ABS(W1).GT.ABS(W2)) THEN
IF(ABS(W1).LE.DE) THEN
GOTO 5
ELSE
FLAG=1
GOTO 4
END IF
END IF
5
GOTO 4
CONTINUE
RETURN
END
Приложение D
SUBROUTINE FFMAX (DT,AL,M,NM,WM)
C
C
1
2
Определение максимального значения
NM-ой производной ИХ формирующей функции
WM=0.
K=0
K=K+1
T=(K-1)*DT
CALL FFSTNM (T,AL,M,NM,FF)
IF (WM.GT.FF) GO TO 2
WM=FF
GO TO 1
CONTINUE
RETURN
END
Приложение E
SUBROUTINE FFSTNM (T,AL,M,NM,FF)
C
C
Вычисление ИХ NM-ой производной
формирующей функции
F0=EXP(-AL*T)
JF=NM-1
47
1
CALL FACT (JF,F1)
FF=AL**JF
S=0.
DO 1 I=1,NM
JI=NM-I
I1=I-1
CALL FACT (JI,F2)
CALL FACT (I1,F3)
B=(-1)**JI*F1/(F2*F3)
MI=M+1-I
CALL FACT (MI,F4)
V=(AL**(M+2-I)/F4)*T**(MI)*F0
S=S+B*V
FF=FF*B*S
RETURN
END
Приложение F
SUBROUTINE FACT (N,F)
C
1
2
Вычисление факториала
F=1
IF (N.EQ.0) GO TO 2
IF (N.EQ.1) GO TO 2
DO 1 I=1,N
F=F*I
RETURN
END
Приложение G
SUBROUTINE FFSPNM (KW,DT,AL,M,NM,NF,OMM)
C
C
Определение максимальной частоты спектра
NM-ой производной ИХ формирующей функции
REAL P(2048),Q(2048)
DO 1 K=1,KW
T=(K-1)*DT
48
1
2
3
4
5
6
CALL FFSTNM (T,AL,M,NM,FF)
P(K)=FF
N1=NF/2+1
CALL N1APFK (KW,P,Q,DT,NF,N1,-1)
SM=0.
DO 2 K=1,N1
P(K)=SQRT(P(K)*P(K)+Q(K)*Q(K))
IF (SM.GT.P(K)) GO TO 2
SM=P(K)
KSM=K
CONTINUE
DOM=6.2831852/(NF*DT)
DO 3 K=1,N1
OM=(K-1)*DOM
P(K)=P(K)/SM
CONTINUE
K=KSM
OM=(K-1)*DOM
IF (P(K).GT.0.05) GO TO 5
KOM=K
OMM=OM
GO TO 6
K=K+1
GO TO 4
CONTINUE
RETURN
END
Приложение H
SUBROUTINE FFIXYX (NK,AL,M,DT,KN,W,NF,N1,AYX)
C
C
1
2
Вычисление производных (0 - N) ИХ
формирующих звеньев
REAL W(KN,NK),AYX(NK,N1),P(2048),Q(2048)
DO 2 K=1,KN
T=(K-1)*DT
DO 1 J=1,NK
CALL FFSTNM (T,AL,M,J,FF)
W(K,J)=FF
WRITE (1,101) T,(W(K,I),I=1,NK)
CONTINUE
49
C
3
4
5
6
15
101
Вычисление АЧХ формирующих функций
PI=3.14
DOM=2.*PI/(NF*DT)
DO 6 I=1,NK
DO 3 K=1,KN
P(K)=W(K,I)
CALL N1APFK (KN,P,Q,DT,NF,N1,-1)
PM=0.
DO 4 K=1,N1
P(K)=SQRT(P(K)**2+Q(K)**2)
IF (PM.LT.P(K)) PM=P(K)
CONTINUE
DO 5 K=1,N1
P(K)=P(K)/PM
AYX(I,K)=P(K)
CONTINUE
DO 15 K=1,N1
OM=(K-1)*DOM
WRITE (2,101) OM,(AYX(I,K),I=1,NK)
FORMAT (1X,F7.3,5(1X,E11.4))
RETURN
END
Приложение I
PROGRAM INTEGRAL
USE IFNLS
Integer CP
C
Объявление переменных
REAL*4 X(1000), XX(1000,5), WF(1000,5), W(1000)
REAL F1(5), F2(5)
CP = NLSSetLocale("Russian", "Russia", 1251)
C
Сравнение значений интегралов F1 и F2
PRINT *, 'Параметры ФФ'
NK=5
PRINT *, 'ALF->'
READ* , ALF
PRINT *, 'MF->'
READ* , MF
50
DT=0.01
C
C
Определение эффективной длительности ИХ
формирующей функции
PRINT*, 'Уровень усечения ИХ DEL->'
READ* , DE
CALL FFTWNM (DT,ALF,MF,NK,DE,KWF,TWF)
PRINT 101,TWF,KWF
C
C
Вычисление ИХ формирующей функции и ее производных
и запись в файл FORT1
DO 2 K=1,KWF
T=(K-1)*DT
DO 1 J=1,NK
CALL FFSTNM (T,ALF,MF,J,FF)
WF(K,J)=FF
WRITE (1,102) T,(WF(K,I),I=1,NK)
CONTINUE
1
2
C
C
3
C
4
5
C
6
Заполнение сигнала и 4-х его производных
и запись в файл FORT3
KON=KWF
DO 3 K=1,KON
T=(K-1)*DT
XX(K,1)=SIN(2*T)
XX(K,2)=2*COS(2*T)
XX(K,3)=-4*SIN(2*T)
XX(K,4)=-8*COS(2*T)
XX(K,5)=16*SIN(2*T)
WRITE (3,102) T,(XX(K,I),I=1,NK)
CONTINUE
Производные сигнала
DO 6 J=2,NK
DO 4 K=1,KON
X(K)=XX(K,J)
DO 5 K=1,KWF
W(K)=WF(K,1)
Интеграл свертки
CALL N1YSWI (KWF,W,KWF,X,Y)
F1(J-1)=Y*DT
CONTINUE
51
C
7
8
9
Производные ФФ
DO 9 J=2,NK
DO 7 K=1,KON
X(K)=XX(K,1)
DO 8 K=1,KWF
W(K)=WF(K,J)
CALL N1YSWI (KWF,W,KWF,X,Y)
F2(J-1)=Y*DT
PRINT 103, F1(J-1),F2(J-1)
CONTINUE
101
FORMAT (2X,'Эффект. длит. ИХ ФФ_NM TW=',F7.3,'
KW=',I4)
102
FORMAT (1X,F7.3,5(1X,E11.4))
103
FORMAT (2X,'F1=',F10.5,' F2=',F10.5)
END PROGRAM INTEGRAL
Приложение J
SUBROUTINE N1YSWI (KW,W,KON,X,Y)
C
Интеграл свертки на конечном интервале времени Tw (KW)
REAL*4 W(KW),X(KON)
1
2
3
4
DO 4 K=1,KON
K1=K+1
Y=0.
IF (K.GT.KW) GO TO 2
DO 1 I=1,K
Y=Y+W(I)*X(K1-I)
GO TO 4
CONTINUE
DO 3 I=1,KW
Y=Y+W(I)*X(K1-I)
CONTINUE
RETURN
END
52
Отзывы:
Авторизуйтесь, чтобы оставить отзыв