сотрудник
Москва, г. Москва и Московская область, Россия
сотрудник
Москва, г. Москва и Московская область, Россия
сотрудник
Москва, г. Москва и Московская область, Россия
сотрудник
Москва, г. Москва и Московская область, Россия
ВАК 05.13.18 Математическое моделирование, численные методы и комплексы программ
УДК 681.5.015 Идентификация, моделирование, параметры и т. п.
ББК 302 Проектирование
Цель исследования заключается в повышении эффективности идентификации параметров динамических систем, имеющих в своем составе упругий элемент с зоной нечувствительности и сухим трением по экспериментальным данным. Задача, решению которой посвящена статья, состоит в разработке алгоритма идентификации, позволяющего определить параметры частотной передаточной функции нелинейной динамической системы по всем экспериментально полученным точкам частотного годографа системы. Методы исследования базируются на теории динамических систем, в частности на гармонической линеаризации. Новизна работы состоит в применении в качестве меры близости для годографов системы и модели суммы квадратов модулей модификации отклонений точек годографа модели системы от соответствующих экспериментально полученных точек годографа системы. Результатом исследования является алгоритм идентификации нелинейной системы, который сводится к решению системы уравнений, линейной относительно неизвестных параметров передаточной функции модели системы. Вычислительный эксперимент показал, что ошибки определения значений параметров системы с зоной нечувствительности упругого элемента и нелинейностью типа «сухое трение» приближенно равны ошибкам измерений экспериментальных отсчетов годографа этой системы. Основной вывод по исследованию состоит в том, что разработанный алгоритм целесообразно применять для решения практических задач идентификации нелинейных систем, в состав которых входят упругий элемент с зоной нечувствительности и сухим трением.
динамическая система, зона нечувствительности упругого элемента, сухое трение, идентификация, частотный годограф
Введение
Первичная разработка, модификация и прототипирование различных динамических устройств связаны с необходимостью проведения большого количества экспериментов. В ходе этих экспериментов характеристики устройств могут быть получены посредством изучения отклика на заранее заданные тестовые воздействия, применяемые к устройствам. Устройства, обладающие рядом заранее известных характеристик, объединяются в соответствующий класс. В настоящей работе рассматривается задача идентификации динамической системы, которой соответствует класс устройств, имеющих в своем составе упругие элементы с зоной нечувствительности и с сухим трением.
Актуальность работы обусловлена необходимостью проведения экспериментальных исследований на этапе разработки, прототипирования и внедрения в производство различных устройств в машиностроении и на транспорте. В настоящее время в этих отраслях в нашей стране возрастает активность коллективов разработчиков во взаимодействии с производителями различного оборудования. Кроме того, в ряде случаев требуется доработка уже эксплуатирующегося оборудования. В частности, в хлопкопрядении к устройствам рассматриваемого класса относятся следящие приводы, которые управляют линейной плотностью выпускаемого полуфабриката (чесальной ленты) [1]. От стабильности линейной плотности продукта, получаемого с чесальных, ленточных и ровничных машин, во многом зависит качество пряжи и, как следствие, качество выпускаемых тканей. Регуляторы линейной плотности на каждом из участков технологического процесса управляют соответствующими вытяжными устройствами. От скорости рабочих органов вытяжного устройства непосредственно зависит линейная плотность получаемого полуфабриката. В состав следящего привода входят упругий элемент с зоной нечувствительности и сухим трением. В производстве проката качество продукции во многом зависит от стабильности ее геометрических размеров. Оборудование прокатных станов содержит приводы клетей, в состав которых также входят устройства рассматриваемого класса [2].
Предполагается, что с реального устройства, описываемого нелинейным дифференциальным уравнением второго порядка, путем подачи гармонического сигнала получены амплитудно-частотная и фазо-частотная характеристики. Решив соответствующие уравнения, можно построить динамическую модель того или иного устройства, которая необходима для его дальнейшего совершенствования. Для решения уравнений необходимо путем обработки экспериментальных данных определить все их коэффициенты. Кроме того, ввиду наличия случайных воздействий в виде помех на реальных производствах, необходимо оценить величину возможных ошибок, обусловленных этими помехами. Такую оценку можно получить на основе статистического моделирования. При синусоидальном входном воздействии на динамическую систему в ней возникнут незатухающие колебания с амплитудой и частотой, которые можно зафиксировать в ходе эксперимента. На основе этих данных возможно применить для идентификации метод гармонической линеаризации. Установившиеся в системе незатухающие колебания будут описываться следующим уравнением:
(1)
В этом уравнении F(x(t)) – зависимость силы, действующей на пружину от координаты x(t). Эта зависимость содержит зону нечувствительности и приведена на рис 1.
Рис. 1. Зависимость силы, действующей на пружину от координаты x(t)
Fig. 1. The dependence of the force acting on the spring on the x(t) coordinate
На рис. 1: k – коэффициент жесткости пружины, b – параметр, определяющий ширину зоны нечувствительности. Зависимость силы трения от координаты 
, определяющая сухое трение, приведена на рис. 2
Рис. 2. Зависимость силы трения от координаты 
,
определяющая сухое трение
Fig. 2. The dependence of the friction force
on the coordinate x ̇(t) determining dry friction
На рис. 2: g - сила сухого трения.
Материалы, модели, эксперименты и методы.
Исходными данными для нахождения значений параметров k, b, a2, g уравнения (1) служат экспериментально полученные амплитудно-частотная и фазо-частотная характеристики системы, представленные в виде годографа. Зависимости силы, действующей на пружину с зоной нечувствительности от координаты x(t) и зависимость силы трения g от координаты , определяющая сухое трение, можно заменить на зависимости, определяемые коэффициентами гармонической линеаризации [3,4]. При такой замене вынужденные колебания системы будут определяться уравнением
, (2)
где, в соответствии с [2]
, (3)
a1(ω) = 4g / (πA(ω)ω). (4)
Здесь A(ω) - амплитуда синусоидальной составляющей колебаний системы на частоте ω. Рассмотрим колебания системы, которые удовлетворяют условию 
. С учетом допущения, что
, 
,
коэффициент гармонической линеаризации e0(ω) для таких колебаний можно вместо (3) описать упрощенным выражением:
a0(ω) = k - 4kb / (πA(ω)) (5)
с учетом обозначения
s = 4kb / π (6)
вместо (5) можем записать
a0(ω) = k - s / A(ω). (7)
С учетом обозначения
v = 4g / π (8)
вместо (4) можем записать
a1(ω) = 4v / (A(ω)ω). (9)
Подставляя в (2) коэффициенты гармонической линеаризации (7) и (9), получим частотную передаточную функцию вида:
. (10)
Параметры k, s, a2, v частотной передаточной функции (10), можно получить, решив соответствующую задачу идентификации [5-8].
Функцию (10) можно переписать в альтернативном виде:
W(jω) = P(ω) + jQ(ω),
где P(ω) и jQ(ω) - вещественная и мнимая части частотной передаточной функции соответственно, определяемые следующими выражениями:
,
(11)
Амплитуду вынужденных колебаний динамической системы для частоты ω можно получить из соотношения A2(ω) = P2(ω) + Q2(ω) или с учетом (11)
.
Здесь R = (k - a2 ω2)2, T = -2(k - a2 ω2)s, U = s2+ v2 – 1.
Частотную передаточную функцию W(jω) можно представить в виде соответствующего годографа на комплексной плоскости.
При экспериментальном оценивании частотных характеристик реальных устройств, особенно в реальных производственных условиях необходимо учитывать возможные случайные воздействия (помехи). Эти воздействия приводят к тому, что точки годографа смещаются случайным образом. Обозначим вещественные и мнимые значения отсчетов экспериментально полученного годографа Wэ(jω) устройства для rexp значений частот ω1, ω2,..., ωrexp
P1 =P(ω1),..., Prexp =P(ωrexp), Q1 = Q(ω1),..., Qrexp = Q(ωrexp) (12)
С учетом того, что вид передаточной функции известен заранее, соответствующий годограф будет принадлежать классу годографов, задаваемых моделью Wм(jω) с учетом (10)
.
или в виде
Wм(jω) = α/( γ + jδ),
где
α = 1, γ = km - sm /A(ω) - e2 ω2, δ = vm /A(ω).
Идентификацию рассматриваемой системы можно свести к решению системы линейных уравнений относительно неизвестных параметров km,, sm, e2, vm [9]
, (13)
.
где
, 
, 
.
Здесь элементы матрицы Φ и элементы вектора 
имеют значения
, 
, 
,
, 
, 
,
,
, 
,
, 
, 
.
Решая систему уравнений (13), получим выражения для параметров km,, sm, e2 в виде формул Крамера
km, = H0 / H, sm, = H1 / H, e2, = H0 / H, где
, 
, 
,
а для параметра vm
С учетом обозначений (6) и (8) получим
bm, = π sm / (4km), gm, = π vm / 4.
Вычислительный эксперимент по оценке ошибок в значениях параметров частотной передаточной функции системы, полученных с применением разработанного алгоритма, проведен для значений k = 1,0; b = 0,1; e2 = 1,0; g = 0,5. Ошибки измерения отсчетов Pi, Qi (12) экспериментально полученного годографа были выбраны случайно в соответствии с равномерным законом плотности распределения вероятностей в интервале [-0,1; 0,1]. Число вычислительных экспериментов rexp = 20. Проводилось rser=200 серий таких экспериментов в частотном диапазоне от 0.2 до 2.0 Гц. В сериях экспериментов вычислялись шибки errk = k - km, errb = b - bm, erra2 = a2 - e2, errg = g - gm, определения параметров k, b, a2 , g и среднеквадратические отклонения sko для случайных величин errk, errb, erra2, errg:
sko_errk =0,0279; sko_errb = 0,0234; sko_errg =,.0150.
Также построены гистограммы для этих величин. На рис. 3-5, приведены примеры гистограмм для случайных величин errk, errb, errg.
|
14 |
Рис. 3. Гистограмма ошибки errk определения параметра k
Fig. 3. Histogram of the error errk for determining the parameter k
Рис. 4. Гистограмма ошибки errb определения параметра b
Fig. 4. Histogram of the error errb of determining the parameter b
Рис. 5. Гистограмма ошибки errg определения параметра g
Fig. 5. Histogram of the error errg for determining the parameter g
Результаты
Проведенный вычислительный эксперимент показал, что ошибки определения значении параметров k, b, g имеют близкий к нормальному закону распределения. На рис. 3-5 приведены соответствующие этим параметрам гистограммы ошибок. Для нормального закона распределения примерно 95% результатов определения значения параметра находится в доверительном интервале ± 2sko и примерно 99% результатов находится в доверительном интервале ±3sko [10]. Оценки погрешностей определения значений параметров k, b, g в вычислительном эксперименте в виде доверительных интервалов приведены в таблице.
Таблица
Доверительные интервалы результатов определения значений параметров
Table
Confidence intervals of the results of determining the values of parameters
|
|
par=k |
par=b |
par=g |
|
sko_errpar (par=k,b,g) |
0.0279 |
0.0234 |
0.0150 |
|
± 2sko_errpar (par=k,b,g) |
±0.0558 |
±0.0468 |
±0.0300 |
|
± 3sko_errpar (par=k,b,g) |
±0.0737 |
±0.0702 |
±0.0450 |
Таким образом показано, что ошибки определения значения параметров k = 1,0; b = 0,1; g = 0,5 в вычислительном эксперименте приближенно соответствуют интервалу ошибки измерений отсчетов годографа [-0,1; 0,1].
Обсуждение/Заключение
Предложен алгоритм идентификации параметров нелинейной динамической системы с зоной нечувствительности упругого элемента и нелинейностью типа «сухое трение» по экспериментально полученным отсчетам частотного годографа системы. Вычислительный эксперимент показал, что ошибки определения значений параметров системы сравнимы с интервалом ошибки измерений экспериментальных отсчетов годографа этой системы. Поэтому алгоритм пригоден для исследования реальных устройств, в состав которых входят рассмотренные нелинейности.
Несмотря на все более широкое внедрение в практику методов искусственного интеллекта, классические методы и алгоритмы анализа динамических систем остаются актуальными и развиваются, многие специалисты продолжают работы в этом направлении [11-13]. Наиболее перспективным по мнению авторов подходом при этом является гибридизация, которая с одной стороны предусматривает использование современных вычислительных возможностей по обработке «больших» экспериментальных данных с помощью нейросетевых [11,12], нечётких [13] и других моделей, а с другой предусматривает основу в виде классических методов и алгоритмов анализа таких систем [14]. Кроме того, существуют гибридные модели на основе теории игр, рассматривающие поведение агентов как динамическую систему в контексте решения социальных дилемм [15]. Предложенный в настоящей статье алгоритм является развитием классических методов анализа динамических систем и может применяться для идентификации реальных устройств в процессе их разработки, модификации и прототипирования, а также в качестве компонента той или иной гибридной модели, предназначенной, например, для решения социальных дилемм со множеством агентов. Предварительная идентификация параметров системы позволит сформировать обоснованные критерии исключения из массива накапливаемых данных заведомо ложных элементов, обусловленных различными сбоями и помехами при их измерении и записи.
1. Севостьянов П.А. и др. Компьютерная модель изменения характеристик волокнистого материала в технологическом процессе. Известия Вузов. Технология текстильной промышленности. 2016. №. 1. С. 361.
2. Мещеряков В.Н., Толчеев В.М. Математическая модель электромеханической системы электропривода прокатной клети стана холодной прокатки. Электротехника: сетевой электронный научный журнал. 2015. Т. 2. №. 2. С. 12-21.
3. Кулябов Д.С., Королькова А.В., Велиева Т.Р. Применение метода гармонической линеаризации к исследованию автоколебательного режима систем с управлением. Discrete and Continuous Models and Applied Computational Science. 2017. Т. 25. №. 3.
4. Попов Е.П., Пальтов И.П. Приближенные методы исследования нелинейных автоматических систем. М.: ГИФМЛ, 1960. 790 с.
5. Леонов Г.А. О методе гармонической линеаризации. Автоматика и телемеханика. 2009. №. 5. С. 65-75.
6. Спиридонов Е.Г., Лещенко М.Н. Разработка алгоритмов оценки параметров динамических систем объектов транспортной инфраструктуры. Транспорт: наука, образование, производство (ТРАНСПОРТ-2021). 2021. С. 211-214.
7. Павлов Ю.Н., Недашковский В.М, Тихомирова Е.А., Шавырин И.Б. Метод гармонической линеаризации в задаче идентификации нелинейных динамических систем. Наука и образование: электронное научно-техническое издание. 2014. № 4. С. 382-394. doi:https://doi.org/10.7463/0414.0704613.
8. Шапкарин А.В., Просандеев А.В., Кулло И.Г. Анализ нелинейных систем автоматического управления методом гармонического баланса в среде MATLAB. Прикаспийский журнал: управление и высокие технологии. 2013. №. 1. С. 077-085.
9. Бойков И.В., Кривулин Н.П. Методы идентификации динамических систем. Программные системы: теория и приложения. 2014. Т. 5. №. 5 (23). С. 79-96.
10. Вентцель Е.С. Теория вероятностей: Учеб. для вузов. 6-е изд. стер. М.: Высш. шк., 1999. 576 с.
11. Tiumentsev Y., Egorchev M. Neural network modeling and identification of dynamical systems. Academic Press, 2019.
12. Quaranta G., Lacarbonara W., Masri S. F. A review on computational intelligence for identification of nonlinear dynamical systems. Nonlinear Dynamics. 2020. Т. 99. №. 2. С. 1709-1761.
13. Škrjanc I., Blažič S., Agamennoni O. Identification of dynamical systems with a robust interval fuzzy model. Automatica. 2005. Т. 41. №. 2. С. 327-332.
14. Семенов А.Д. и др. Комбинированная нейросетевая система регулирования линейной плотности ленты. Известия высших учебных заведений. Технология текстильной промышленности. 2021. №. 2. С. 109-112.
15. Akiyama E., Kaneko K. Dynamical systems game theory II: A new approach to the problem of the social dilemma. Physica D: Nonlinear Phenomena. 2002. Т. 167. №. 1-2. С. 36-71.
