Афанасьев А.А. Математические основы теории систем управления
Подождите немного. Документ загружается.
n
nn
x
h
xx
λ
=
−
+
~~
1
или
.
~
)(
~
)1(
~
0
1
1
xhxhx
n
nn
+
+
+=+=
λ
λ
Отсюда видно, что условие
01
xxx
nn
≤≤≤
+
K
будет выполнено, если
11
≤+
h
λ
. (4.24)
Поэтому справедливо 02111
≤≤−⇒≤+≤−
hh
λ
λ
. Или
20
≤−≤
h
λ
. (4.25)
Таким образом, явная схема Эйлера условно устойчива, т.е. устойчива,
если шаг интегрирования удовлетворяет условию
λ
2
≤
h
или
⎟
⎟
⎠
⎞
⎜
⎜
⎝
⎛
=≤
λ
1
2 TTh , (4.26)
где
Т- постоянная времени в решении уравнения (4.23).
При численном интегрировании системы уравнений состояния тестовое
уравнение имеет вид
[]
.
xA
d
t
dx
=
(4.27)
Как известно, эту систему можно расщепить на n независимых уравнений
(см.п.2.16):
[]
0
1
0
,,,2,1, xSznkz
d
t
dz
kk
k
−
===
K
λ
(4.28)
путем замены переменных
[ ]
,zSx
=
где [S]- матрица правых собственных векторов матрицы [A] .
Поскольку собственные числа систем (4.27) и (4.28) совпадают, то для
оценки устойчивости разностных схем можно рассматривать отдельные
уравнения (4.28). Если собственные значения
K
λ
вещественные, то
неравенство (4.26), ограничивающее шаг интегрирования, запишется так :
131
max
h
λ
2
≤
. (4.29)
При наличии комплексных собственных значений
)(
kkkk
j
ω
α
λ
λ
+=
неравенство (4.24) получит вид (после возведения в квадрат левых и правых
частей его)
1)()1(
22
≤++
hh
kk
ω
α
. (4.30)
Это неравенство задает круг единичного радиуса с центром в точке (рис.4.2):
,1 0
=−=
h
k
h
k
ω
α
.
Множество комплексных значений ,
удовлетворяющих условию асимптотической
устойчивости решения разностной схемы для
тестового уравнения, называют областью
устойчивости метода.
λh
Обратим внимание, что при чисто мнимом
)0(
=
kk
α
λ
область устойчивости (4.30) для явного
метода Эйлера вырождается в нуль.
jh
K
ω
α
K
h
-2 -1 0
Рис . 4.2
Из (4.30) можем получить неравенство, ограничивающее шаг
интегрирования при комплексном значении
k
λ
02)(121
2222222
≤++⇒≤+++
hhhh
kkkkkk
α
ω
α
ω
α
α
.
Отсюда
2
22
Re22
k
k
kk
k
h
λ
λ
ωα
α
−=
+
−≤
.
Это неравенство и условие (4.29) указывает, что при интегрировании (4.27) с
большим по модулю собственным значением шаг интегрирования должен быть
выбран достаточно малым. При этом попытка увеличить шаг более его
максимальной оценки приводит к резкому возрастанию погрешности-” взрыву”
погрешности.
2) Неявная схема Эйлера
При применении этого метода для решения тестового уравнения (4.23)
получим
132
1
1
~
~~
+
+
=
−
n
nn
x
h
xx
λ
или
h
x
x
n
n
λ
−
=
+
1
~
~
1
.
Поскольку при вещественном 0
<
λ
1
1
1
<
−
h
λ
при любых 0
>−
h
λ
, то схема безусловно устойчива.
Для комплексных собственных значений
kkk
j
ω
α
λ
+=
условие
устойчивости
1
1
1
<
−
h
k
λ
(4.31)
может быть последовательно преобразовано к видам :
1)()1(11
22
>+−⇒>−
hhh
kkk
ω
α
λ
.
Последнее неравенство задает внешность единичного круга с центром в точке
0;1
==
hh
kk
(рис. 4.3) .
α
ω
jh
k
ω
α
k
h
0
1
2
-1
1
Рис. 4.3
Областью устойчивости неявного
метода является вся плоскость за
исключением упомянутого выше
единичного круга. В этом круге 0
>
h
k
α
,
т.е.
0
>
k
α
и, следовательно, тестовое
уравнение неустойчиво.
Видим, что при интегрировании
абсолютно устойчивых уравнений
состояния условие устойчивости (4.31) не
накладывает каких-либо ограничений на
шаг дискретизации. Выбор шага в этом
случае должен осуществляться только по
соображениям точности вычислений.
Обратим внимание,что решения разностного уравнения неявного метода
Эйлера оказываются устойчивыми и в правой полуплоскости (за исключением
единичного круга), где решение исходного дифференциального уравнения
неустойчиво.
Следовательно, использование данного метода для интегрирования
неустойчивых уравнений состояния дает неверный результат, не
соответствующий характеру истинного решения.
133
3) Метод трапеций
При решении тестового уравнения (4.23) разностным методом трапеций
(4.9) будем иметь
)
~~
(
2
~~
11 nnnn
xx
h
xx
++=
++
λ
.
Этот результат получаем из формулы (4.9) после подстановки:
x
x
f
~
)
~
(
λ
=
.
После преобразования разностная схема примет вид
nn
x
h
h
x
~
21
21
~
1
λ
λ
−
+
=
+
Из этого равенства также следует, что
0
1
1
~
21
21
~
x
h
h
x
n
n
+
+
⎟
⎟
⎠
⎞
⎜
⎜
⎝
⎛
−
+
=
λ
λ
. (4.32)
Для обеспечения устойчивости разностной схемы необходимо выполнение
условия
1
21
21
<
−
+
=
h
h
q
λ
λ
(4.33)
Или
()( ) ()( )
2222
221221 hhhh
kkkk
ωαωα
+−<++ 0<⇒−<⇒
kkk
hh
α
α
α
. (4.34)
Областью устойчивости этого метода является
вся левая полуплоскость
(рис.4.4) .
jh
k
ω
0
α
k
h
Рис. 4.4
При чисто мнимом значении )(
ω
λ
λ
j=
устойчивы как тестовое уравнение, так и
рассматриваемый разностный метод, так как в этом
случае, как видно из (4.33) , q =1. Это обстоятельство
позволяет использовать метод трапеций при
интегрировании уравнения (4.27), собственные
значения матрицы [A] которого содержат пару
мнимых сопряженных величин. Решение в этом случае
134
имеет составляющую незатухающего гармонического колебания.
4) Схема Рунге-Кутта второго порядка.
Решим тестовое уравнение (4.23) с помощью разностной схемы (4.19).
Выполняя в этой схеме подстановку
x
x
t
f
~
)
~
,(
λ
=
или
⎥
⎦
⎤
⎢
⎣
⎡
++===
−
••
+
)
~
,(
2
~
,
2
),(
2
1
nnnnnnn
nn
xtf
h
x
h
tfxtfk
h
xx
.
Опустим в правой части этой формулы время
⎟
⎠
⎞
⎜
⎝
⎛
+⇒ )
~
(
2
~
xf
h
xf
n
.
Далее вместо функции оставим только её аргумент, как требует вид тестового
уравнения (4.23)
⎟
⎠
⎞
⎜
⎝
⎛
+=⇒
nn
x
h
xr
~
2
~
λλ
получим
nnnnnn
xqx
h
xhxrhxx
~~
2
~~~~
1
=
⎟
⎠
⎞
⎜
⎝
⎛
++=+=
+
λλ
,
где
()
2
1
2
h
hq
λ
λ
++= .
Схема устойчива, если
( )
1
2
1
2
<++=
h
hq
λ
λ
. (4.35)
Для действительных отрицательных
λ
требование (4.35) выполняется
при
2<h
λ
. (4.36)
Это условие совпадает с ограничением на шаг для явного метода Эйлера.
Внутренняя область, ограниченная кривой 2 на рис 4.5, является
областью устойчивости рассматриваемого метода применительно к
комплексным значениям
λ
(номер кривой на рис. 4.5 совпадает с порядком
точности метода).
135
5) Схема Рунге-Кутта четвертого порядка
Решение тестового уравнения (4.23) применительно к разностной схеме
(4.21) будет иметь вид
nn
x
hhh
hx
~
2462
1
~
443322
1
⎟
⎠
⎞
⎜
⎝
⎛
++++=
+
λλλ
λ
,
а соответствующее ограничение на шаг задается неравенством
1
2462
1
443322
<++++
hhh
h
λλλ
λ
. (4.37)
Для вещественного 0>
λ
неравенство
(4.37) дает оценку
7853,2≈< Ch
λ
,
а для чисто мнимого )(
ω
λ
λ
j= имеем
22<h
λ
.
Кривая 4 на комплексной
плоскости рис.4.5 ограничивает область
устойчивости этого метода.
6) Многошаговые методы Адамса
.
а) Явные методы Адамса
Для явных методов Адамса
ограничения на шаг более
значительные, чем в методах Рунге-Кутта [14]:
4
2
1
-3
-2
-1
0
1
2
3
-1
-2
-3
jh
k
ω
α
k
h
Рис. 4.5
r
Ch <
λ
.
Метод первого порядка точности (метод Эйлера) имеет ;
2
1
=
C
второго - 1
2
=C ;
третьего - 556,0
3
≈C ;
четвертого -
3,0
4
=
C
.
136
С дальнейшим увеличением порядка точности величина продолжает
убывать.
r
C
Для сравения отметим, что для действительных отрицательных
λ
в
методе Рунга-Кутта четвертого порядка шаг мог быть выбран почти в десять
раз большим, чем в аналогичном методе Адамса.
Области устойчивости (верхняя симметричная половина) явных методов
Адамса для комплексных значений
λ
( номер кривой совпадает с порядком
точности метода) приведены на рис. 4.6.
Анализ этих областей и областей на рис 4.5. вынуждает забраковать
традиционные явные методы
Рунге-Кутта и Адамса
применительно к решению
жестких систем, так как при
r
Ch <
λ
и больших
λ
использование малых значений h
требуют значительных затрат на
вычисления.
Для жестких уравнений
необходимы такие методы, у
которых область устойчивости
при решении тестового уравнения
(4.23) содержала бы всю левую
полуплоскость плоскости h
λ
, т.е.
при 0Re <
λ
решение соответствующего разностного уравнения было бы
асимптотически устойчиво при любом положительном h.
hj
k
ω
0.8
1
2 0,4
3 4
2,0 -1,6 1,2 -0.8 -0,4 0
h
k
α
Рис. 4. 6
Эти методы получили название А-устойчивых, и для них оказались
справедливыми следующие утверждения:
1) никакие явные методы не могут быть А-устойчивыми;
2) не существует А-устойчивого неявного линейного многошагового метода с
порядками точности 3 и выше.
б) Неявные методы Адамса
Неявные методы Адамса первого (Эйлера) и второго (Адамса с m=1)
порядка точности уже были рассмотрены выше. Это неявные методы
соответственно Эйлера и трапеций.
Областью устойчивости первого метода будет вся плоскость h
λ
за
исключением единичного круга с центром (1,0) и, следовательно, неявный
метод Эйлера А-устойчив.
Для второго метода (неявный метод Адамса с m=1 и p=2) областью
устойчивости является вся левая (отрицательная) полуплоскость )0Re(
<h
λ
, т.е.
метод трапеций также является А-устойчивым. Однако при 1−<<h
λ
137
быстрозатухающему решению (4.23)
0
~~
xex
hn
n
λ
= соответствует решение
метода трапеций (4.32)
0
~
2
1
~
x
h
x
n
⎜
⎝
⎛
+=
λ
2
1
h
nn −
⎟
⎠
⎞
⎜
⎝
⎛
−
⎟
⎠
⎞
λ
,
которое отнюдь не является быстроубывающим, так как при
−∞→h
λ
имеем
0
~~
xx
n
−→
.
Таким образом, требование А-устойчивости не является достаточным для
гарантии “хороших” свойств методов, и для метода трапеции шаг
интегрирования не может быть выбран слишком большим из-за невозможности
получения необходимой точности решения.
Неявный метод Эйлера свободен от этого недостатка, так как
−∞→→ hприx
n
λ
0
~
.
Неявные методы Адамса более высоких порядков точности, согласно
высказанному выше утверждению 2 , не могут быть А-устойчивыми.
Для 0<h
λ
справедливы ограничения :
r
C
метод третьего порядка (4.16) (с числом шагов m=2) имеет
6
3
=< Ch
λ
;
метод четвертого порядка (4.17) (m=3) –
3
4
=< Ch
λ
(как и у метода Рунге-Кутта четвертого порядка).
Для того чтобы быстроубывающей компоненте решения
дифференциального уравнения
(
− → − ∞λ
k
отвечала быстроубывающая
компонента решения разностного уравнения, используются такие А-устойчивые
методы, у которых дробно- рациональная функция
)( hq
λ
в решении тестового
уравнения
nn
xhqx
~
)(
~
1
λ
=
+
имеет степень полинома в знаменателе выше, чем степень в числителе, и
0)(lim =
∞−→
hq
h
λ
λ
.
138
Такие методы получили название L-устойчивых. К их числу относятся неявные
(не классические) методы Рунге- Кутта [14] и неявный метод Эйлера.
У него
h
hq
λ
λ
−
=
1
1
)(.
Тот факт, что явные классические методы (Эйлера, Адамса, Рунге- Кутта)
не А-устойчивы, свидетельствует о том, что их использование для
интегрирования некоторых уравнений состояния приводит к большим
вычислительным трудностям. На существование подобных уравнений, трудно
поддающихся интегрированию явными классическими методами, впервые
обратили внимание в 1952 г. и назвали их жесткими системами
дифференциальных уравнений.
В настоящее время существует специальная теория жестких уравнений и
методов их решения [14].
Заметим, что применительно к уравнениям электрических цепей
жесткость является скорее правилом, чем исключением.
Вопросы для самопроверки
1.
Как выглядят тестовые скалярные и векторно-матричные уравнения
состояния, которые используются для оценки устойчивости численных
методов при варьировании начальных условий?
2.
Назовите неравенства, ограничивающие шаг интегрирования явным
методом Эйлера при вещественных и комплексных значениях
собственных чисел матрицы тестового уравнения состояния.
3.
Каковы области устойчивости явного и неявного методов Эйлера?
4.
Назовите область устойчивости метода трапеций.
5.
Какова область устойчивости
A
- устойчивых методов?
6.
Могут ли быть
A
- устойчивыми явные методы?
7.
Какие методы называются
L
- устойчивыми?
4.4. Жесткость уравнений состояния
4.4.1
.
Понятие жесткости
Рассмотрим дифференциальное уравнение первого порядка с малым
параметром
μ
при производной
0),,( >=
μμ
xtf
d
t
dx
(4.38)
139
для случая, когда вырожденное уравнение, соответствующее (4.38),
0)
~
,( =
x
t
f
(4.39)
имеет единственное решение
)()(
~
t
G
t
x
= (4.40)
и в окрестности этого решения величина
xf
∂
∂
отрицательна. Последнее
условие является достаточным для устойчивости решения
)()(
~
t
G
t
x
=
.
x (t
0
,x
0
)
=G(t)
0
t
Рис. 4.7
Для начальной (стартовой) точки (t
0
. ), расположенной выше кривой
0
x
)()(
~
t
G
t
x
= , имеем для производной
0)()(
),(1),()())((
00
0
0
<−=−≅==
−
=
xxxx
x
xtfxtf
dt
tdx
dt
xtxd
xx
μ
λ
∂
∂
μμ
,
где
0
),(1
0
<=
=xx
x
xtf
∂
∂
μ
λ
.
Если обозначить
0
)()( xtxty −= ,
то справедливо
t
eyty
μ
λ
0
)( = .
140
~
x