close

Вход

Забыли?

вход по аккаунту

?

Аналог метода Хилла в задачах устойчивости периодических решений нелинейных дискретно-континуальных систем.

код для вставкиСкачать
ПРОБЛЕМЫ ЕСТЕСТВЕННЫХ НАУК
УДК 531.36
Д.К. Андрейченко
АНАЛОГ МЕТОДА ХИЛЛА В ЗАДАЧАХ УСТОЙЧИВОСТИ ПЕРИОДИЧЕСКИХ
РЕШЕНИЙ НЕЛИНЕЙНЫХ ДИСКРЕТНО-КОНТИНУАЛЬНЫХ СИСТЕМ
Предложен аналог метода Хилла исследования устойчивости
периодических решений применительно к дискретно-континуальным
механическим системам, движение дискретных элементов которых
описывается
нелинейными
обыкновенными
дифференциальными
уравнениями с запаздывающим аргументом.
D.K. Andreichenko
THE HILL METHOD ANALOG IN THE STABILITY PROBLEMS
OF PERIODIC SOLUTIONS OF NONLINEAR DISCRETE-CONTINUAL SYSTEMS
The analog of the stability Hill method research of periodic solutions
with reference to discrete - continual mechanical systems is offered here,
driving of which discrete elements are described by the nonlinear ordinary
differential equations with lagging argument.
α
Введение
..
L
yc
C
Исследование устойчивости в малом периодических
a
решений
нелинейных
систем
обыкновенных
O L
дифференциальных уравнений, в соответствии с теорией Lc a
0
Флоке, сводится к вычислению матрицы монодромии и
N
0
y
исследованию ее собственных значений [1, с.187].
Применительно к системам обыкновенных дифференциальных
уравнений с запаздывающим аргументом либо к дискретноz
континуальным
системам
(ДКС),
описываемым
y(z,t)
совокупностями обыкновенных дифференциальных уравнений
O1 l
и связанными с ними граничными условиями и условиями
N
1
связи уравнениями в частных производных, численное
y(t)
1
1
построение
матрицы
монодромии
становится
L1
проблематичным. В настоящей статье, на примере
рассмотренной ранее в [2] нелинейной газореактивной
α1
системы стабилизации спутника с упругим стержнем и
z
закрепленным на его конце телом, показана принципиальная
Рис. 1
возможность развития метода Хилла [3, с. 269] к исследованию
устойчивости периодических решений ДКС, движение
континуальных элементов которых описывается линейными уравнениями в частных
производных с не зависящими явно от времени коэффициентами, а движение дискретных
элементов – нелинейными обыкновенными дифференциальными уравнениями с
запаздывающим аргументом.
1. Уравнения движения дискретно-континуальной системы
Уравнения движения показанной на рис. 1 ДКС (совершающего плоское движение
абсолютно жесткого спутника с моментом инерции Jc и массой mc под действием
возмущающего момента L, с жестко закрепленным на расстоянии a от центра массы
спутника прямолинейным однородным упругим стержнем, несущим жестко закрепленное
на противоположном конце абсолютно жесткое тело с массой m1 и моментом инерции J1)
в безразмерной форме имеют вид
t
&& = L − f ( w(t − τ)) + L0 − aN 0 , w(t ) = b1α
& (t ) + b2 α(t ) + b3 ∫ α( ξ)d ξ ,
J cα
(1.1)
&& + α
&&1 ) = L1 , m1[ &&
&& ] = N1 ,
mc &&
yc = N 0 , J 1 ( α
yc + &&
y1 − (1 + a )α
(1.2)
&& = −( y ′′′′( z , t ) + γy& ′′′′( z, t )), ( )′ = ∂ / ∂z ,
&&
y ( z, t ) + &&
y c − ( a + z )α
(1.3)
0
z = 0 : y (0, t ) = y ′(0, t ) = 0;
z = 1: y (1, t ) = y1 (t ), y′(1, t ) = −α1 (t ) ,
(1.4)
L0 = − y ′′(0, t ) − γy& ′′(0, t ), N 0 = − y ′′′(0, t ) − γy& ′′′(0, t )
L1 = y ′′(1, t ) + γy& ′′(1, t ), N1 = y ′′′(1, t ) + γy& ′′′(1, t )
(1.5)
t ≤ 0 : α(t ) = α& (t ) = yc (t ) = y& c (t ) = y ( z, t ) = y& ( z, t ) = 0 .
(1.6)
Здесь f(w) – некоторая функция своего аргумента, описывающая нелинейность типа
насыщения, b1, b2, b3 – коэффициенты обратных связей ( )′ = ∂ ( ) / ∂z , ( &) = ∂ ( ) / ∂t .
Примем далее
f ( w) = th( w) .
Линеаризация в окрестности состояния равновесия приводит уравнение (1.1) к виду
(уравнения (1.2)-(1.6) при этом своего вида не изменяют)
&& = L − b1α& (t − τ) − b2α(t − τ) − b3 ∫
J cα
t −τ
0
α( ξ)d ξ + L0 − aN 0 , j = 1,2,3 .
(1.1.а)
Области устойчивости линейной ДКС (1.1.а), (1.2)-(1.6) в пространстве параметров
обратных связей b1, b2, b3, приведенные в [5], подробно исследовались на основе
предложенных в [4] теорем об устойчивом квазимногочлене и об устойчивых,
неустойчивых и асимптотически устойчивых квазирациональных дробях, а также
предложенного в [5] варианта метода D-разбиений. Следует отметить, что наличие малого
запаздывания τ резко ограничивало области устойчивости за счет дестабилизации
управляемой деформируемой конструкции по высшим формам колебаний. В случае, когда
параметры обратных связей выходят за пределы областей устойчивости, исходные
уравнения (1.1)-(1.6) при L(t)=0 (с точностью до несущественного выбора начального
момента времени) допускают 2π/ω-периодические решения
α = α(0) (t ), α1 = α1(0) (t ), yC = yC(0) (t ), y1 = y1(0) (t ), y ( z, t ) = y (0) ( z, t ) ,
(1.7)
которые описывают возможные периодические автоколебательные режимы с частотой ω,
исследованные численно в [2] на основе метода гармонического баланса в высших
приближениях. Во всех случаях период возможных автоколебаний (даже по высшим
формам) превосходил характерное время запаздывания:
2π / ω > τ .
(1.8)
При численном моделировании развития начальных возмущений [2] выполнялась
дискретизация уравнения (1.1) на основе проекционного метода Бубнова-Галеркина с
использованием в качестве координатных базисных функций ортогональных полиномов
Чебышева 1-го рода, и далее полученная система обыкновенных дифференциальных
уравнений с запаздывающим аргументом интегрировалась численно. При этом было
выяснено, что на фазовой плоскости ( α, α& ) фазовые траектории достаточно быстро
попадали в окрестность одного из возможных предельных циклов, однако вопрос о
строгом исследовании устойчивости предельных циклов по Ляпунову остался открытым.
2. Исследование устойчивости предельных циклов
С точностью до выбора несущественного начального момента времени t=0,
линейные уравнения движения ДКС, возмущенной в малой окрестности возможного
предельного цикла, имеют вид (уравнения (1.2)-(1.5) своего вида не изменяют)
&& = − f (0) (t − τ)w(t − τ) + L0 − aN 0 , w(t ) = b1α& (t ) + b2 α(t ) + b3v(t ), v&(t ) = α(t )
J cα
(
)
f (0) (t ) = f ′ b1α& (0) (t ) + b2 α(0) (t ) + b3 ∫ α(0) (t )dt , f ′( z ) = 1/ ch 2 ( z )
(2.1)
Обозначим
Y (t ) = {v (t ), α(t ), α& (t ), α1 (t ), α& 1 (t ), yc (t ), y& c (t ), y1 (t ), y&1 (t ), y ( z, t ), y& ( z, t )} .
Аналогично [1, с. 298], поскольку исходные уравнения (1.1), (1.2)-(1.5) не зависят
от времени, с точностью до несущественного выбора начального момента времени при
L(t ) ≡ 0 уравнения возмущенного движения (2.1), (1.2)-(1.5) по форме записи аналогичны
уравнениям относительно набора величин
&& (t ), α& (t ), α
&& (t ), y& (t ), &&
Y& (t ) = {v&(t ) = α(t ), α& (t ), α
y (t ), y& (t ), &&
y (t ), y& ( z, t ), &&
y ( z, t )} ,
1
1
c
c
1
1
полученным дифференцированием по времени (1.1), (1.2)-(1.5). Следовательно, 2π/ωпериодические функции
v (t ) = α(0) (t ), α = α& (0) (t ), α1 = α& 1(0) (t ), yC = y& C(0) (t ), y1 = y&1(0) (t ), y ( z, t ) = y& (0) ( z, t )
(2.2)
заведомо будут решением уравнения возмущенного движения.
Модельные уравнения (2.1), (1.2)-(1.5) задают некоторый линейный оператор
Lt : Y |[ −τ ,0] → Y |[ t −τ,t ]
(2.3)
(строго говоря, задание «предыстории» на отрезке [ −τ,0] требуется лишь для функций
v (t ) = ∫ α(t )dt , α(t ), α& (t ) , для остальных функций требуется задание начальных условий).
Далее, поскольку правые части линейных обыкновенных дифференциальных уравнений
(2.1), (1.2) являются достаточно гладкими функциями времени, а уравнение (1.3)
учитывает малое, но ненулевое вязкое трение в материале конструкции согласно модели
Фойгта, относительно негладким «предысториям» Y |[ −τ,0] при t > τ будут соответствовать
более гладкие решения Y |[ t −τ,t ] (рассматриваемые как функции пространственной
переменной и времени). Следовательно, оператор Lt при t > τ будет вполне непрерывен в
смысле [6].
Так как коэффициент линейного уравнения (2.1) при слагаемом с запаздывающим
аргументом является 2 π / ω -периодической функцией времени, при выполнении условия
(1.7)
Lt +T = Lt L T = L T Lt , T = 2π / ω
(2.4)
Lt + nT = Lt L Tn = L Tn Lt , T = 2π / ω, n = 1, 2,3,...
(2.5)
и соответственно
Следует отметить, что условие (1.8) не является достаточно жестким
ограничением; так, если оно не выполняется, следует выбрать период T>τ, кратный
периоду автоколебаний, и соответственно уменьшить частоту ω.
Из (2.5) непосредственно следует, что, в соответствии с теорией Флоке, возможные
частные решения уравнений (2.1), (1.2)-(1.5) будут представлять собой произведение
показательной функции времени t на 2 π / ω -периодические функции времени
iωnt
α(t ) = eλt ∑ n =−∞ αn eiωnt , α1 (t ) = eλt ∑ n =−∞ α(1)
, yC (t ) = eλt ∑ n =−∞ yn( C )eiωnt
ν e
∞
∞
∞
y1 (t ) = eλt ∑ n =−∞ yn(1)eiωnt , y ( z, t ) = eλt ∑ n =−∞ yn ( z )eiωnt
∞
∞
(2.6)
Далее, легко убедиться, что совокупность значений функций (2.6) при t ∈ [ −τ,0]
представляет собой «собственный элемент» линейного оператора L T , T = 2π / ω , а
величина eλT = e 2 πλ / ω является соответствующим ему собственным значением. Поскольку
оператор L T вполне непрерывен, его спектр ограничен, состоит из не более чем счетного
множества собственных значений и имеет единственную возможную точку сгущения
eλT = 0 [6, с. 189], чему соответствует Re λ → −∞ . Таким образом, имеется лишь конечное
число решений, которым соответствуют показатели λ, лежащие в правой комплексной
полуплоскости Re λ>0, и обеспечивающие неустойчивость предельных циклов. Из (2.2)
следует наличие среди (2.6) по крайней мере одного чисто периодического решения,
которому соответствует показатель λ=0. Следует отметить, что добавление к показателю λ
величины iω приводит к сдвигу на 1 индекса Фурье-коэффициентов в (2.6), т.е.
фактически подобные показатели неразличимы.
Далее полагаем
f (0) (t ) = ∑ n =−∞ f n(0)eiωnt .
∞
(2.7)
Выражения (2.2) позволяют удовлетворить уравнениям движения (1.2)-(1.5) (в силу
(c)
(1)
их линейности), при этом уравнения относительно αn , α(1)
n , yn , yn , yn ( z ) аналогичны
уравнениям (2.2)-(2.5) из [2] при замене λ → λ + inω и приводятся к следующему
результату
(1)
2
(c)
ϕν1 (λ + inω)αn + ϕν 2 (λ + inω)α(1)
n + ϕ ν 3 ( λ + inω) yn + ( λ + inω) ϕ ν 4 ( λ + inω) yn = 0
n = 0, ±1, ±2,..., ν = 2,3, 4
(2.8)
Подстановка (2.6), (2.7) в (2.1) и перемножение рядов Фурье приводит к
следующим линейным уравнениям (вообще говоря, бесконечным)
(1)
2
(c)
ϕ11 (λ + iωn )αn + ϕ12 (λ + iωn )α(1)
n + ϕ13 (λ + iωn ) yn + ( λ + i ωn ) ϕ14 ( λ + iωn ) yn =
= −(λ + iωn )e −τ ( λ+iωn ) ∑ j ≠n f n(0)
n = 0, ±1, ±2,...
− j [ b1 ( λ + iωj ) + b2 + b3 /( λ + iωj )]α j ,
Здесь
ϕ11 (λ ) = ( J c + aμ21 − μ11 )λ 3 + f 0(0) (b1λ 2 + b2λ + b3 )e −τλ ,
а остальные величины полностью аналогичны [2]
ϕ14 (λ ) = λ( aμ 24 − μ14 ) , ϕ12 (λ ) = (1 + γλ )λ( aμ 22 − μ12 ) , ϕ13 (λ ) = (1 + γλ )λ(aμ 23 − μ13 ) ,
ϕ21 (λ ) = −μ 21λ 2 , ϕ22 (λ ) = −(1 + γλ )μ 22 , ϕ23 (λ ) = −(1 + γλ )μ 23 , ϕ24 (λ ) = mc λ 2 − μ 24 ,
ϕ31 (λ ) = ( J 1 − μ31 )λ 2 , ϕ32 (λ ) = J 1λ 2 − (1 + γλ )μ 32 , ϕ33 (λ ) = −(1 + γλ )μ33 , ϕ34 (λ ) = −μ34 ,
ϕ41 = −[m1 (1 + a ) + μ 41 ]λ 2 , ϕ42 (λ ) = −(1 + γλ )μ 42 , ϕ43 (λ ) = m1λ 2 − (1 + γλ )μ 43 , ϕ44 = m1λ 2 − μ 44 ,
μ νj = μ νj [k (λ )], ν, j = 1, 2,3, 4 , k (λ ) = eiπ / 4 λ (1 + γλ ) −1/ 4 ,
(2.9)
(2)
μ ν1 ( k ) = aμ(1)
ν1 ( k ) + μ ν1 ( k ), ν = 1, 2,3, 4 ,
(1)
(2)
μ11
( k ) = k −2 (ch k − cos k − sh k sin k ) / σ( k ) , μ11
( k ) = k −3[k (ch k − cos k ) − sin k (ch k − 1) +
+ sh k (cos k − 1)]/ σ(k ) , μ12 ( k ) = k (sh k − sin k ) / σ( k ) , μ13 (k ) = k 2 (ch k − cos k ) / σ(k ) ,
(1)
−1
μ14 ( k ) = −μ11
( k ) , μ(1)
21 ( k ) = k [sin k (ch k − 1) + sh k (cos k − 1)] / σ( k ) ,
(2.10)
−2
μ(2)
21 ( k ) = k [ch k − cos k − k (sh k + sin k ) + sh k sin k ] / σ( k ) , μ 22 ( k ) = −μ13 ( k ) ,
μ 23 ( k ) = −k 3 (sh k + sin k ) / σ(k ) , μ 24 ( k ) = −μ (1)
21 ( k ) ,
(1)
(2)
−3
μ(1)
31 ( k ) = −μ11 ( k ) , μ 31 ( k ) = k [(1 − ch k )sin k + (cos k + k sin k − 1)sh k ] / σ( k ) ,
(1)
μ32 ( k ) = k (ch k sin k − sh k cos k ) / σ( k ) , μ33 ( k ) = k 2 sh k sin k / σ(k ) , μ34 ( k ) = μ11
(k ) ,
−2
μ(2)
41 = k (cos k − ch k + k ch k sin k + ( k cos k − sin k )sh k ) / σ( k ) ,
(1)
3
(1)
μ(1)
41 ( k ) = μ 21 ( k ) , μ 42 ( k ) = μ 33 ( k ) , μ 43 ( k ) = k (ch k sin k + sh k cos k ) / σ( k ) , μ 44 ( k ) = −μ 21 ( k ) ,
σ( k ) = ch k cos k − 1 .
Функции μ νj ( k ), ν, j = 1, 2,3, 4 аналитичны в окрестности точки k=0 и при конечных
γ>0 не имеют особенностей в правой половине и на мнимой оси комплексной плоскости
(λ). Легко проверить
lim λ
λ→∞
−β νj
μ νj [k (λ )] = bνj = const , Reλ ≥ 0, ν, j = 1, 2,3, 4 ,
(2.11)
b11 = a γ1/ 2 , β11 = −1/ 2 , b12 = b13 = 0 , β12 = β13 = 0 , b14 = −γ −1/ 2 , β14 = 1/ 2 , b21 = (1 + i )a γ1/ 4 / 2 ,
β21 = −1/ 4 , b22 = b23 = 0 , β22 = β23 = 0 , b24 = 2 γ −3/ 4 , β24 = 3/ 4 , b31 = −i 2 γ 3/ 4 , β31 = −3/ 4 ,
b32 = − 2 γ −1/ 4 , β32 = 1/ 4 , b33 = b34 = γ −1/ 2 , β33 = β34 = 1/ 2 , b41 = ( a + 1) 2 γ1/ 4 , β41 = −1/ 4 ,
b42 = γ −1/ 2 , β42 = 1/ 2 , b43 = b44 = 2 γ −3/ 4 , β43 = β44 = 3/ 4 ,
μ νj [k (λ )] = μ νj [k (λ )] .
(2.12)
3. Аналог метода Хилла
Вводя обозначение vν = α ν /(λ + iων) , удобно представить
однородную систему линейных уравнений (2.8), (2.9) в форме
бесконечную
vν + Φ (λ + iων)e − τ ( λ+iων ) ∑ j ≠ν f ν−(0)j ⎡⎣β1 (λ + iωj ) 2 + β 2 (λ + iωj ) + β3 ⎤⎦ v j = 0, ν = 0, ±1, ±2,... (3.1)
Здесь
Φ (λ ) = Q (λ ) / D(λ ) , D(λ ) = det{ϕνj (λ )}, ν, j = 1, 2,3, 4 ; Q (λ ) = det{ϕνj (λ )}, ν, j = 2,3, 4
(3.2)
Φ (λ ) ,
D (λ ) ,
Q (λ ) –
соответственно,
аналоги
передаточной
функции,
характеристического и возмущающего квазимногочленов ДКС [2], [4]. В соответствии с
(2.10)-(2.12), функции D(λ ) и Q (λ ) аналитичны при Re λ ≥ 0 , и
Imλ
D(λ ) = D(λ ), Q (λ ) = Q (λ ),
λ
lim λ −7 D (λ ) = Const , Re λ ≥ 0
λ→∞
(3.3)
Рассуждения, аналогичные п. 2, позволяют
установить, что характеристический квазимногочлен
Reλ
Рис. 2
D(λ) имеет не более чем конечное число нулей при Re λ>0. Заметим, что аналитичность
функций D(λ) и Q(λ), а также конечность числа нулей функции D(λ) в правой
полуплоскости позволяют достаточно быстро на основании принципа аргумента провести
проверку наличия у них совпадающих нулей в правой комплексной полуплоскости. Для
этого достаточно сначала проверить наличие нулей D(λ) и Q(λ) в некоторых вертикальных
полубесконечных полосах, а затем, если таковые будут выявлены, проверить наличие
нулей D(λ ) и Q (λ ) в прямоугольниках некоторой конечной высоты в пределах
полуполосы (рис. 2).
Из приведенного в [4] доказательства теоремы об устойчивом квазимногочлене
следует, что, в соответствии с (3.3), число нулей (с учетом их кратности)
характеристического квазимногочлена D(λ) в правой полуплоскости Re λ>0 определяется
величиной
1
(3.4)
P=
7 − 2 Δ arg D (iσ) .
0≤σ≤∞
2π
Будем далее предполагать, что квазимногочлены D(λ) и Q(λ) не имеют
совпадающих нулей при Re λ>0; в таком случае величина (3.4) определяет суммарную
кратность полюсов передаточной функции Ф(λ) в правой полуплоскости.
Определитель бесконечной системы линейных уравнений (3.1), являющийся
аналогом определителя Хилла, вычисляется следующим образом:
(
)
D (λ ) = lim Dn (λ ), Dn (λ ) = det{Dνj (λ )}, − n ≤ ν ≤ n, − n ≤ j ≤ n ,
(3.5)
Dνν (λ ) = 1, Dνj (λ ) = f ν−(0)j e−τ ( λ+iων )Φ(λ + iων) ⎡⎣b1 (λ + iωj )2 + b2 (λ + iωj ) + b3 ⎤⎦ , ν ≠ j .
(3.6)
n →∞
Используя неравенство Адамара [7, с. 35], можно показать, что в некоторой
конечной области правой комплексной полуплоскости, последовательность (3.5) сходится
равномерно, и, следовательно, функция D (λ ) обладает свойствами, типичными для
определителя Хилла:
D (λ ) = D (λ ) ,
(3.7)
D (λ + iω) = D (λ ) ,
(3.8)
D (λ ) → 1 при Re λ → ∞ .
(3.9)
Сверх того, из равномерной сходимости последовательности (3.5) и свойства
периодичности (3.8) следует, что функция D (λ ) будет аналитична при Re λ ≥ 0 , за
исключением счетного числа полюсов, получающихся сдвигом полюсов передаточной
функции Ф(λ) (т.е. нулей характеристического квазимногочлена D(λ)) на величину inω ,
n = 0, ±1, ±2,...
С вычислительной точки зрения, последовательность (3.5) может сходиться
достаточно медленно, однако сходимость соответствующего ей ряда
D (λ ) = 1 + ( D1 (λ ) − 1) + ( D2 (λ ) − D1 (λ )) + ( D3 (λ ) − D2 (λ )) + ...
значительно улучшается преобразованием в непрерывную дробь [8].
Можно показать, что, если в правой части бесконечной однородной системы
уравнений (3.1) заменить нули некоторой последовательностью, для которой сходится
сумма квадратов модулей элементов, то решение можно найти по формулам,
аналогичным формулам Крамера, и оно также будет представлять последовательность,
для которой будет сходиться сумма квадратов модулей элементов. Следовательно, для
того, чтобы бесконечная однородная система линейных уравнений (3.1) допускала
нетривиальное решение относительно Фурье-коэффициентов
vν , ν = 0, ±1, ±2,... ,
необходимо потребовать выполнения условия, аналогичного уравнению Хилла
D (λ ) = 0 .
(3.10)
Аналогично теореме Андронова-Витта можно полагать, что предельный цикл будет
орбитально асимптотически устойчив [1, с. 309], если, с точностью до несущественного
сдвига на кратную iω величину, все возможные показатели λ будут иметь отрицательные
мнимые части, кроме одного показателя λ=0 кратности 1.
Применительно к обыкновенным дифференциальным уравнениям, когда
отсутствует запаздывание аргумента, а передаточная функция Ф(λ) представляет собой
рациональную дробь, расположение и кратность полюсов которой известны, имеется
возможность выделения особенностей определителя Хилла D (λ ) и его дальнейшего
упрощения на основе теоремы Лиувилля. Однако в случае ДКС сложность поведения
передаточной функции Φ (λ ) в левой комплексной полуплоскости исключает подобную
возможность.
Imλ
l1
Imλ
Imλ
λ
λ
λ
2iω
l2
ω
0
Reλ
а
iω
ω
0 Reλ
)
б
0
Reλ
)
в
Рис. 3
Из свойства периодичности (3.8) следует, что, при исследовании устойчивости
предельного цикла достаточно исследовать наличие либо отсутствие корней определителя
Хилла D (λ ) в некоторой полубесконечной полосе (рис. 3, a) в правой комплексной
полуплоскости, параллельной действительной оси и имеющей ширину ω, причем там
может быть лишь конечное число корней. Далее, из свойства периодичности следует, что
в рассматриваемую полуполосу попадет лишь столько различных полюсов определителя
Хилла, сколько имеется в правой полуплоскости различных полюсов передаточной
функции Ф(λ), если только полюса не отличаются друг от друга на величину, кратную iω.
С другой стороны, из (3.6) непосредственно следует, что, если λ=λ0 и λ=λ1 являются
различными полюсами передаточной функции Ф(λ), отличающимися на величину,
кратную iω , то λ=λ0 будет полюсом определителя Хилла, причем его кратность будет
равна суммарной кратности двух указанных полюсов передаточной функции.
Следовательно, суммарная кратность полюсов определителя Хилла в рассматриваемой
полубесконечной полосе априорно известна, совпадает с суммарной кратностью полюсов
передаточной функции Ф(λ) в правой полуплоскости и определяется величиной (3.4).
Однако в подобном случае число нулей в рассматриваемой полуполосе легко
устанавливается на основе принципа аргумента.
Значение λ=0 априорно является корнем определителя Хилла D (λ ) , причем
устойчивому периодическому предельному циклу соответствует корень λ=0 кратности 1.
Следовательно, с учетом (3.8), (3.9), суммарная кратность полюсов P и «опасных» нулей N
определителя Хилла в полубесконечной полосе (пути обхода контуров l1 и l 2 показаны
на рис. 3, a и б соответственно) оценивается выражениями
P=
1 ⎛
⎞
7π − 2 Δ arg D (λ ) ⎟ ,
⎜
2π ⎝
λ∈l 2
⎠
(3.11)
(
)
1
(3.12)
2πP − Δ arg D (λ ) − 1 .
λ∈l1
2π
Полуполоса на рис. 3, a должна выбираться так, чтобы ее верхняя и, следовательно,
нижняя полубесконечные границы не проходили через полюса определителя Хилла D (λ ) ,
т.е. через нули функций D(λ + inω), n = 0, ±1, ±2,... . Такой выбор расположения верхней
границы заведомо возможен, т.к. в правой комплексной полуплоскости в любой
параллельной действительной оси полубесконечной полосе ширины ω функция D (λ )
имеет не более чем конечное число полюсов, и может быть достигнут, например,
перемещением верхней границы по вертикали (с некоторым шагом) от точки iω / 2 до
точки iω . В свою очередь, проверка отсутствия нулей функций D (λ + inω), n = 0, ±1, ±2,...
в узкой окрестности верхней границы проверяется посредством анализа изменения
аргумента указанных функций при обходе контура, показанного на рис. 3, в.
В случае, когда входящая в (2.1) функция f (0) (t ) мало отличается от константы,
т.е. коэффициенты соответствующих модельных уравнений практически не зависят от
времени, для коэффициентов Фурье-разложения (2.7) выполняется условие f v( 0 ) << f 0( 0 ) ,
N=
ν = ±1, ±2, ±3,... В тех областях правой комплексной полуплоскости, где отсутствуют
полюса передаточных функций Φ (λ + inω), n = 0, ±1, ±2,... , внедиагональные элементы
определителя Хилла (3.5) будут пренебрежимо малы. Следовательно, некоторые из
возможных нулей определителя Хилла будут содержаться в достаточно малых
окрестностях полюсов передаточных функций Φ (λ + inω), n = 0, ±1, ±2,... , т.е. в малых
окрестностях корней (нулей) характеристического квазимногочлена D(λ ) , сдвинутых на
величину, кратную iω .
В качестве примера рассмотрим спутник с безразмерными
Таблица 1
параметрами a=0,05, mc=3, Jc=0,5, m1=1, J1=0,3, γ=0,08, τ=0,11,
b1=0,5, b2=16,5, b3=5. В данном случае параметры обратных связей
Re λ
Im λ
не принадлежат области устойчивости, а приведенные в табл. 1
-0,3030
0
собственные частоты колебаний позволяют сделать вывод о
-0,0314
1,281
дестабилизации
3-й
формы
колебаний
деформируемой
-0,2979
4,600
конструкции. Система допускает три предельных цикла с
0,2734
5,989
частотами автоколебаний ω=1,25, 4,39, 5,81, показанные на
рис. 4, a штриховой линией, пунктиром и сплошной линией соответственно. Фазовые
траектории (показаны на рис. 4, б пунктиром) на плоскости ( α, α& ) , соответствующие
процессу развития малых начальных возмущений, достаточно быстро подходят к
предельному циклу с частотой ω=5,81 (показан сплошной линией). Приращение
аргумента характеристического определителя D (λ ) при обходе контура l 2 , показанного
на рис. 3, б, составляет 7π / 2 , а его годограф в специальном масштабе
u + iv = D(λ )ln(1+ | D(λ ) |) / | D(λ ) | приведен на рис. 4, в. Следовательно, согласно (3.11), в
правой полуплоскости корни характеристического определителя D(λ ) отсутствуют.
Далее, приращение аргумента определителя Хилла D (λ ) при обходе контура l1 ,
показанного на рис. 3, a, равно 2π , а его годограф показан на рис. 4, г. Согласно (3.12),
определитель Хилла не имеет «опасных» корней в полосе, показанной на рис. 3, a, и
предельный цикл с частотой ω=5,81 устойчив по отношению к малым возмущениям.
Ближайшим к мнимой оси корнем определителя Хилла D (λ ) , кроме λ = 0 , оказывается
λ = -0,0368 ± 1,41i .
.
α
.
α
ω=4,39
0,15
0,2 ω=1,25
ω=5,81
0
0
-0,15
-0,2
-0,4
-0,14 -0,07
0
a)
а
0,07
α
v
-0,3
-0,05 -0,025
0
0,025 α
b) б
ImD
0
0,8
-40
0
-0,8
-80
-100
-10
ω=5,81
-5
0
u
-1,6
c)
в
0
1
2
3 ReD
d) г
Рис. 4
Таблица 2
Два других предельных цикла с частотами автоколебаний
ω=1,25 и ω=4,39 оказываются неустойчивыми. Согласно (3.12),
Re λ
Im λ
для предельного цикла с частотой ω=1,25 в показанной на рис. 3, a
-0,9948
0
полубесконечной полосе определитель Хилла D (λ ) кроме корня
-0,00076
1,003
λ=0 имеет еще 3 корня ( λ = 0,0561 и λ = 0,1463 ± 0,6002i ). Для
-0,8905
13,69
предельного цикла с частотой ω=4,39 определитель Хилла
-5,550
34,05
согласно (3.12) имеет 1 «опасный» корень (λ=0,4024).
6,861
50,55
В качестве другого примера приведем рассмотренный
ранее в [2,5] спутник с параметрами a=0,05, mc=34,75, Jc=0,07442, m1=3, J1=0,007, γ=0,01,
τ=0,01, b1=0,7, b2=200, b3=150. В данном случае параметры обратных связей не
принадлежат области устойчивости, а приведенные в табл. 2 собственные частоты
колебаний позволяют сделать вывод о дестабилизации 4-й формы колебаний
деформируемой конструкции. Возникает предположение о том, что развитие малых
возмущений может привести к возникновению автоколебаний.
Формально, система допускает три предельных цикла с частотами автоколебаний
ω=6,787, 11,46, 20,12, показанные на рис. 5, a штриховой линией, пунктиром и сплошной
линией соответственно. В начальные моменты времени 0 ≤ t ≤ 4,5 фазовые траектории
(показаны на рис. 5, б пунктиром) на плоскости ( α, α& ) , соответствующие процессу
развития малых начальных возмущений, достаточно быстро подходят к предельному
циклу с частотой ω=20,12 (показан сплошной линией). Далее, как показано на рис. 5, в, в
моменты времени 4,5≤t≤180 (чему соответствует примерно 560 предполагаемых периодов
автоколебаний), фазовые траектории заполняют некоторую относительно узкую
окрестность предельного цикла, и становится неясно, является ли предельный цикл
устойчивым. Приращение аргумента характеристического определителя D (λ ) при обходе
контура l 2 , показанного на рис. 3, б, составляет –π/2, а его годограф приведен на рис. 5, г.
.
α ω=6,787
.
α ω=20,12
3
0,7
ω=20,12
0
0
ω=11,46
-0,7
-3
-6
-0,8
-0,4
0
а a)
.
α
0,4
α
-1,4
-0,06 -0,03
0,03
α
v
0
0,7
0
-40
-0,7
-80
-1,4
-0,07 -0,035
0
в c)
-100
-24
0,035 α
ImD
.
α
15
3
0
0
-15
-3
-30
-10
0
b)
б
0
10 20
д e)
30 ReD
-6
-40
-12
0
u
0
α
d)г
-20
f)е
Рис. 5
Следовательно, согласно (3.11), D(λ ) имеет две пары комплексно сопряженных
корней в правой полуплоскости ( λ = 0,05583 ± 0,9162i
и λ = 0,1251 ± 10,19i
соответственно). Далее, приращение аргумента определителя Хилла D (λ ) при обходе
контура l1 , показанного на рис. 3, a, равно нулю, а его годограф показан на рис. 5, д.
Согласно (3.12), определитель Хилла имеет в полосе, показанной на рис. 3, a, две пары
комплексно сопряженных корней ( λ = 0,04654 ± 0,9373i и λ = 0,03037 ± 9,136i
соответственно), которые и обеспечивают постепенное отталкивание фазовых траекторий
от неустойчивого предельного цикла. Далее, с возрастанием времени t, фазовые
траектории постепенно уходят на бесконечность (рис. 5, е). Следует отметить, что корни
определителя Хилла D (λ ) в рассматриваемом случае оказались достаточно близкими к
корням характеристического квазимногочлена D (λ ) , что в принципе позволило бы
установить слабую неустойчивость возможного автоколебательного режима при помощи
анализа частотного годографа характеристического квазимногочлена D(λ ) .
Два других предельных цикла с частотами автоколебаний ω=6,787 и ω=11,46 также
оказываются неустойчивыми.
Заключение
Следует отметить, что, после несущественных изменений, предложенный метод
будет применим и к исследованию устойчивости вынужденных периодических колебаний
ДКС с нелинейными дискретными элементами. Как показано в настоящей работе,
информация о расположении областей устойчивости управляемых деформируемых
конструкций в пространстве параметров обратных связей имеет существенное значение,
так как малоамплитудные высокочастотные автоколебания, вызванные неудачным
выбором параметров проектируемой конструкции, могут, вообще говоря, оказаться
неустойчивыми, что приведет к дальнейшему (теоретически – неограниченному)
нарастанию начальных возмущений. В связи с этим далее представляет интерес
применение моделей упругих элементов конструкций, более точно описывающих их
поведение в высокочастотной области (например, с использованием моделей
тонкостенных конструкций типа Тимошенко), а также уточнение математических моделей
запаздывающих звеньев в системе управления. Также представляет интерес развитие
строгих методов исследования устойчивости периодических режимов и для ДКС с
нелинейными континуальными элементами. Однако решение данных задач выходит за
рамки настоящей работы.
Работа выполнена при финансовой поддержке Совета по грантам Президента
Российской Федерации для государственной поддержки молодых российских ученых и по
государственной поддержке ведущих научных школ (грант № МД-2328.2005-1).
ЛИТЕРАТУРА
1. Демидович Б.П. Лекции по математической теории устойчивости /
Б.П. Демидович. М.: Изд-во МГУ, 1998. 480 с.
2. Моделирование
влияния
запаздывающего
аргумента
в
нелинейной
газореактивной системе стабилизации спутника с упругим стержнем и закрепленным на
его конце телом / К.П. Андрейченко, Д.К. Андрейченко, В.С. Шорин, С.Г. Наумов /
Вестник СГТУ. 2005. № 3. С. 17-27.
3. Уиттекер Э.Т. Курс современного анализа. Ч.2. Трансцендентные функции /
Э.Т. Уиттекер, Дж.Н. Ватсон. М.: Физматгиз, 1963. 516 с.
4. Андрейченко Д.К. К теории комбинированных динамических систем / Д.К.
Андрейченко, К.П. Андрейченко // Известия РАН. Теория и системы управления. 2000. №
3. С. 54-69.
5. Андрейченко Д.К. К теории стабилизации спутников с упругими стержнями /
Д.К. Андрейченко, К.П. Андрейченко // Известия РАН. Теория и системы управления.
2004. № 6. С. 149-162.
6. Люстерник Л.А. Краткий курс функционального анализа / Л.А. Люстерник,
В.И. Соболев. М.: Высшая школа, 1982. 271 с.
7. Корн Г. Справочник по математике для научных работников и инженеров /
Г. Корн, Т. Корн. М.: Наука, 1978. 832 c.
8. Андрейченко Д.К. Эффективный алгоритм численного обращения интегрального
преобразования Лапласа / Д.К. Андрейченко // Журнал вычислительной математики и
математической физики. 2000. Т. 40. № 7. С. 1030-1044.
Андрейченко Дмитрий Константинович –
доктор физико-математических наук,
профессор кафедры «Механика деформируемого твердого тела»
Саратовского государственного технического университета
1/--страниц
Пожаловаться на содержимое документа