ЖУРНАЛ ВЫЧИСЛИТЕЛЬНОЙ МАТЕМАТИКИ И МАТЕМАТИЧЕСКОЙ ФИЗИКИ, 2007, том 47, < 1, с. 34-38
УДК 519.624
УСТОЙЧИВОСТЬ МАТРИЧНОЙ КОНЕЧНО-РАЗНОСТНОЙ СХЕМЫ ПРИ РЕШЕНИИ ЗАДАЧИ О КОЛЕБАНИЯХ УПРУГИХ СИСТЕМ
© 2007 г. В. Д. Жесткая
(681031 Комсомольск-на-Амуре, пр. Ленина, 27, КнАГТУ) e-mail: office@knastu.ru Поступила в редакцию 25.05.2004 г. Переработанный вариант 21.07.2006 г.
Исследуется устойчивость матричной конечно-разностной схемы, полученной с использованием центральных разностей при численном решении дифференциального уравнения колебаний упругой системы. Библ. 3.
Ключевые слова: метод конечных разностей, устойчивость, дифференциальные уравнения колебаний упругой среды.
Работа посвящена исследованию устойчивости матричной конечно-разностной схемы, полученной с использованием центральных разностей:
- 2ur
M
*г + 1
r-1
*г + 1 '
r - 1
h1
2 h
u0 u-
+ Kur
f 0, f _
r = 0, 1,..., N-1,
(1)
■*-1 ~ '-ъ
для дифференциального уравнения
Ми + Си + Ки = Р, (2)
где М, С и К - вещественные симметричные квадратные матрицы порядка п, а и и Р - векторы (матрицы-столбцы), являющиеся функциями времени л В частности, такой вид имеют уравнения метода конечных элементов для колебаний упругих систем (см. [1]).
Доказательство устойчивости схемы (1) при п = 1 дано в [2], там же предложено получить его для п = 2. Для случая любого п доказательство устойчивости схемы (1) автору неизвестно. Для доказательства первое уравнение в (1) приведем к виду
M + 2 CI ur +1 = - (h K-2 M) ur - IM-2 C) ur-1
Prh2.
(3)
С целью избежать громоздких формул проследим ход рассуждений на примере случая п = 3. Найдем сначала матрицу, обратную к матрице, стоящей в левой части уравнения (3):
M+2C1 =
m11+2 c11 m 12 + 2 c12 m13 + 2 c13
m21 + 2c21 m22 + 2C22 m23 + ^c23
m31 + 2C31 m32 + ^Сз2 m33 + ^C33
/ - - - - \
^11 Д12 ^13 h s11 S12 S13 h2 '11 '12 '13
Ц21 ^22 ^23 + h 2 s21 s22 s23 + — 4 t21 '22 '23
I ^31 ^32 S31 S32 s33 /31 '32 1333 У
а = аег| м + 2 с |,
Цц =
Ш22 Ш23 т32 Ш33
Д12 = -
Ш12 Ш13
Ш32 Ш33
? ... ,
Д33 =
Ш11 Ш12 Ш21 Ш22
(4)
*11 = С22 Ш23 Ш22 С23 , ^12 С12 Ш13 Ш12 С13
+
С32 Ш33 Ш32 С33 С32 Ш33 Ш32 С33
3 II С11 Ш12 + Ш11 С12
С21 Ш22 Ш21 С22
'п =
С22 С23 , '12 = - С12 С13
С32 С33 С32 С33
tзз
С11 С12
С21 С22
Введем следующие обозначения: если А и В - некоторые матрицы, то [AiBjB|¡\ - матрица, первый столбец которой есть ьй столбец матрицы А, второй и третий столбцы - соответственно, j-й и к-й столбцы матрицы В. Аналогично определяются матрицы [BiAjBk\, [BiBjAk\ и определители \AiBjBk |, \BAjBк |, \BBjA,к |.
С учетом этих обозначений получим
а = |м| + 2 (| С1 М2М3\ + \М1 С2М3\ + \Ы1Ы2С3\) +
+ -4 (М1С2 С3 + С М2 С3 + |С1 С2М3) + -8-1 с = 1М + аН + 0( Н), а = 2 (| СМ2М3| + |М1 С2 М3\ + М1М2С3|). Для приведения исследуемой схемы к каноническому виду положим
(5)
Уг
Тогда из (1) и (3) получим
Уг +1 = КиУг + Н Рг
Уо
f о
f-1
оператор шага RН имеет вид
Рг = Н
М + 2С| Р,
-1
(6)
Я,
-| м + 2 сJ (и к -2М) -^м + сJ ^м - с
Е 0
-1
(7)
2
3
и
г
и
г - 1
0
Найдем произведения матриц, входящие в (7):
-2| М| + кЬ11 + к2 р11 + О (к3) М21+ к'Р2!+ О(к3 ) кЬ31 + к2 р31 + О (к3)
М + 2 С | (к К -2 М) =
кЬ 12 + к2Р12 + О(к3) - 2| М| + кЬ22 + к 2 Р22 + О (к3)
кЬ32 + к2Р32 + О(к3 )
кЬ13 + к2 р13 + О (к3) кЬ 23 + к2 Р23 + О (к3) - 2 [ М ] + кЬ33 + к2 Р33 + О (к3)
Г М + к С
.-1
М - '2 С1 =
Здесь
|М| + кd11 + к2 д11 + О (к3) кd12 + к2 д12 + О (к3)
23
кd13 + к д13 + О(к )
кd21+ кд21 + О(к ) |М| + к^22 + к д22 + О(к ) к^23 + к д23 + О(к )
кd31 + к д31 + О (к )
кd32 + к2 д32 + О (к2) |М| + кd33 + к2д33 + О (к3)
Ьц = -2а + |С 1 Ы2Ы3\, Ь; = |СММ3, ] = 2, 3, Ь22 = -2а + МСМ3, Ь2; = С;М3, j = 1, 3, Ь33 = -2а + ММ2С31, Ь3; = ММ2С}\, j = 1, 2,
Ру = X М-гАу - ^ X ^т,-,-, к^ - элементы матрицы К,
г = 1 г = 1
d 11 = а - |С 1 М2М3\, dl] = -IСМ2М3, j = 2, 3,
d22 = а - М1С2М3, d2J = -|МСМ3|, ] = 1, 3,
d33 = а - ММ2С3, dз; = -|М^М2С;|, ] = 1, 2,
33
д>з = 4 Xг
игЩ X
Заметим,что
Ьц + 4,. + а = 0 (. = 1 Подставив (8) и (9) в (7), получим
Як
г = 1 Г = 1
, 3) и Ьч + ^ = 0
Л11 Л 12 Л 13 Л 14 Л 15 Л16
Л21 Л22 Л23 Л24 Л25 Л26
Л31 Л32 Л33 Л34 Л35 Л36
1 0 0 0 0 0
0 1 0 0 0 0
0 0 1 0 0 0
где
Ли = ±(21М| - кЬи + О(к2)), . = 1, 2, 3, Ач =-Х(кЬ1} + О(к2)), . Ф;, j = 1, 2, 3,
(9)
г+3 = (IМ| + мхх + о(И2)), 7 = 1,2,3,
, + о(и2)), 7 Ф], 7, ] = 1, 2, 3.
¿7, ] + 3 д^"-I]
Нормы и(И) и/(И) в выражении Ьиы(К) = /(И) (см. [2]) определим равенствами
||и(и)||ии = тах(|и^|), где и() есть 7-й элемент иг, 7 = 1, 2, 3,
II /-(и)|| _ II/ К =
Рг f 0
f-1
= тах(|/0°|,|, 7 = 1, 2, 3,
а норму в пространстве УИ - в виде
к = тах(|1> Ы, Ы,
|у 1- |у2- УЗ| 1У3- Ув\
При сделанном выборе норм соблюдаются условия теоремы об устойчивости (см. [2, гл. 5]), согласно которой при их выполнении для устойчивости достаточно, чтобы нормы степеней оператора шага были равномерно ограничены по И, т.е. чтобы выполнялось условие
Щу < С3, к = 1, 2,..., N. Если норма в пространстве У задана формулой
||у||у = тах(|У11, ),
то
где
= И»У1 у
(10)
1 0 0 0 0 0
0 1 0 0 0 0
0 0 1 0 0 0
8 = -1 1
и-1 0 0 -И 10 0
0 и-1 0 0 - и-1 0
0 0 И1 0 0 -И-1.
Норма ||КИ||У может быть, с учетом (10), выражена формулой
11КИ1Уи = 8Ки8 У.
Вычислим матрицу SRИS 1. Произведя умножение, получим
¿11 + ¿14 А 12 + ¿15 ¿13 +¿16 ^¿14 ^¿15 -И¿16
¿21 + ¿24 ¿22 + ¿25 ¿23 + ¿26 ^¿24 -И¿25 -И¿26
¿31 + ¿34 ¿32 + ¿35 ¿33 + ¿ 36 -и¿34 И¿ 35 ^¿36
= и (АП-1+А14) 1 ( ^12 + ¿15 ) И(¿13+ ¿ 16) - ¿14 ¿15 - ¿16
И(¿21 + ¿24) 1 ( ^22-1 + ¿25 ) И (¿23+ ¿ 26 ) - ¿24 - ¿25 - ¿26
И ( ¿31 + ¿34) 1 ( ^32 + ¿35 ) 1(¿33- 1+ ¿36 ) - ¿34 - ¿35 - ¿36
у
и
После подстановки значений найдем, что
SRb S —
ви ... вг
в 61 ... в
66
где элементы, стоящие на главной диагонали, имеют вид
12
В. = 1 + кун + — О(к ), причем v11 = ч22 = ч33 = 0,
а остальные элементы - вид
Bl} — I (hvy + О (h1)).
Параметры V., V, входящие в формулы В., В у, есть некоторые числовые величины, зависящие от значений элементов матриц М, С, К. К примеру, ч51 = -\М1К1М31, ч44 = -|М1С2М31, ч56 = -|М1С3М31.
При выбранной в У норме векторов норма матрицы SRhS-1 задается формулой
4-1
SRh S — max
X1В1'I' XlB2iI' XlB3iI' XlB4iI' XВ5г'I' XlB|
Vi — 1
i — 1
i — 1
i—1
i — 1
6i i—1
Структура элементов В ¡у (¡, у = 1, ..., 6) дает возможность выбрать некоторую не зависящую от к постоянную С таким образом, чтобы выполнялось неравенство
1 + Ch > max
X |B1il '.'XIB
\i — 1
6i
i — 1 /
Тогда
R
ш,,
— l|SRhS Y < 1 + Ch,
откуда следует неравенство для степеней оператора шага:
Nik <||Rh||Yh <( 1 + Ch)N < eC, k — 1' 2' ...' N,
что и доказывает устойчивость исследуемой схемы (см. [2], [3]).
Все закономерности, обнаруженные в случае n = 3, очевидным образом переносятся и на n > 3. Так, например, при любом n > 3 выражение для a в формуле (5) примет вид
a — 2(|C1M2...Mn\ + M1C2...Mn\ + ... + M1...Mn_1Cn\),
коэффициент b11 окажется равным -2a + |CjM2.. .Mn |, и т.д. Это позволяет утверждать, что устойчивость схемы (2) при n > 3 доказывается аналогично случаю n = 3, и сделать вывод о том, что (2) устойчива при любом n.
СПИСОК ЛИТЕРАТУРЫ
1. Постное В.А. Численные методы расчета судовых конструкций. Л.: Судостроение, 1977.
2. Годунов С.К., Рябенький В С. Разностные схемы. М.: Наука, 1977.
3. Марчук ГИ. Методы вычислительной математики. М.: Наука, 1980.
Для дальнейшего прочтения статьи необходимо приобрести полный текст. Статьи высылаются в формате PDF на указанную при оплате почту. Время доставки составляет менее 10 минут. Стоимость одной статьи — 150 рублей.