О выборе краевых условий для дифференциально-алгебраических уравнений второго порядка

Ряд природных и технических явлений можно описать взаимосвязанными обыкновенными дифференциальными уравнениями (ОДУ) первого и второго порядков и конечномерными (трансцендентными) уравнениями [8; 10; 14]. Объединяя эти уравнения в виде системы, мы получим систему ОДУ второго порядка с тождественно вырожденной матрицей перед второй производной. Такие системы принято называть дифференциальноалгебраическими уравнениями (ДАУ) второго порядка.

В настоящее время бурно развивается качественная теория и численные методы решения ДАУ первого порядка. Это связано с тем, что такие уравнения имеют принципиальные отличия от систем ОДУ и, кроме того, имеют многочисленные приложения [8-10; 12; 14].

С ДАУ порядка выше первого, как правило, поступают следующим образом: вводят новую вектор-функцию, размерность которой равна kn ( n — размерность исходной задачи, k — порядок ДАУ), и переписывают исходную задачу в виде ДАУ первого порядка [13].

Описанный подход значительно усложняет как качественное исследование, так и создание численных алгоритмов для полученной задачи. Этот факт и инициировал исследования, приведенные ниже, которые являются продолжением работ [3; 4; 7; 11].

1 Постановка задачи

Рассмотрим систему

A(t)х"(t)+ B(t)x\)) + C(t)x(t) = f (t), t e [0,1],            (1)

где A ( t ), B ( t ) и C ( t ) — ( n x n ) матрицы, f ( t ) — заданная, а x ( t ) — неизвестная n -мерные вектор-функции.

В дальнейшем изложении предполагается, что det A (t) = 0.                                   (2)

Системы (1) с условием (2) будем называть дифференциальноалгебраическими уравнениями второго порядка (ДАУ2).

Предполагается, что для системы (1) заданы краевые условия

а о x (0)+ во x '(0)= Y о а1 x (1) + в1 x '(1) = Y1 ’

где а0, в0 — (mо х n) -матрицы, а1, в1 — (m1 х n) -матрицы. В частности, начальную задачу для (1) можно переписать в виде (3) с матрицами

^а о =

Л Z7 Л o’)

V о 7

во =

: ■ ।

V on 7

а , и в 1 — нулевые матрицы, а O_n — единичная

матрица, о — нулевая матрица размерности ( n х n ) .

В отличие от ОДУ второго порядка для выделения единственного ре шения ДАУ2 начальные или краевые условия можно не задавать.

Например, ДАУ2 вида

"ai(t) 0Yxi(t)) .

V о оA x 2 ( t ) 7⁺

[bot)

'

о Y x, (t)) „ I I +

0 7 Vx 2(t )7

Го 1 Y xi(t)

Г о

V1 о 7V x 2 (t )FV ^ )/

t е [ о,1 ] (4)

при достаточно гладких функциях a 1 ( t ), b 1 ( t ) , ^ ( t ) имеет единственное

решение: x1 (t) = ^(t), x2 (t) = -a1 (t)^"(t) - b1 (t)^'(t)

Краевые условия (3) должны также быть заданы корректно, т. е. быть согласованными с правой частью (1).

Приведем критерий корректности задания краевых условий (3). Нам понадобится известное определение.

Определение 1. (см, например, [2; 5]). Матрица A называется полу-обратной к матрице A , если она удовлетворяет матричному уравнению

AA^A = A .

Из этого определения следует, что

VA = (O - AA-)A = о.

Умножаем систему (1) на матрицу V в точках t = о и t = 1.

В силу (5) получим:

' V ( о) С (о) x ( о) + V ( о) B (о ) x'(о) = V (о) f (о), / V(1)С(1)x(1) + V(1)B(1)x'( 1) = V(1) f (1).

Вводя обозначения у _о = ( x T ( о ), x' T ( о )), у1 = ( x ^T ( 1 ), x ' T ( 1 )), V ( 0 ) С ( 0 ) = L _о, V ( о ) B ( 0 ) = M о , V ( о ) f ( 0 ) = N о , V ( 1 ) С ( 1 ) = Д, V ( 1 ) B ( 1 ) = M 1 , V ( 1 ) f ( 1 ) = N 1 ,

краевые условия (3) и условия (6) перепишем в виде:

Д ^а о ^в о) у о = Y о _^{( а} 1 ^в 1 ) У 1 = Y 1 , ⁽ L о M о ) y о = N о \ ( L 1 M 1 ) у 1 = N 1 .

Здесь и всюду в дальнейшем изложении обозначение ( а\в )

означает

расширенную матрицу. Таким образом, если задача (1), (3) имеет решение

(возможно, не единственное), то множество решений системы линейных алгебраических уравнений (СЛАУ) (7) является и множеством решений СЛАУ (8). Обратное неверно.

Для доказательства этого факта достаточно рассмотреть обыкновенное линейное дифференциальное уравнение

x "(t) + B (t)x '(t)+ C (t ) = f (t), где B(t) и C(t) — функции. В этом случае A- — 1, V(0) — V(1) — L0 = — M0 — L1 — M1 — N0 — N1 — 0, то есть решением системы (8) является любой вектор из R2.

2 Методы решения краевой задачи

Опишем алгоритм численного решения краевой задачи (1), (3). Стандартным образом введем на отрезке интегрирования сетку t, — ih , i — 0,1,..., N , h — -1 и обозначим A i — A ( t i ), B i — B ( t i ), C i — C ( t i ) , f ^— f ^{( t}, )> x “ ^{x ( t}, )

Классическая трехточечная дискретизация [1; 6] системы (1) дает

D i x i + 1 + E i x i + F i x i - 1 — G i ,       i — 1,2,..., N - 1,

где ( n x n ) -матрицы D i , E i и F i зависят от формулы, применяемой для аппроксимации (1).

Краевые условия (3) после двухточечной аппроксимации примут вид: J D _O x + E 0 x 0 — G _o,                            (10 a )

[ ^E N x N ^{+ F} N x N - 1 ^{— G} N ,                        ⁽¹⁰ b )

где D ₀ и E ₀ — ( m ₀ x n ) -матрицы, E _N и F _N — ( m ₁ x n ) -матрицы.

Условий (10 a ) , (10 b ) недостаточно для того, чтобы применить метод матричной прогонки для решения СЛАУ (9), (10 a ) , (10 b ) , так как эта система является недоопределенной (содержит ( m ₀ + m ₁ + n ■ ( N - 2 )) уравнений, а неизвестных n ■ N ). Поэтому мы предлагаем поступить следующим образом. Умножим (1) на матрицу V ( t ) — E - A ( t ) A ^-( t ) и вычислим полученное выражение в точках t — 0 и t — 1.

В силу формулы (4) имеем:

V (0)C (0)x (0) + V (0)B (0)x ‘(0) — V (0)f (0)

и

V (1)C (1)x (1)+ V (1)B (1)x ’(1)— V (1)f (1)

Аппроксимируя эти условия, например простейшей формулой численного дифференцирования, получим:

( 11 a ) ( 11 b )

S 0 x 1 + T 0 x 0 - H 0 ,

SNxN + TNxN-1 — HN, где So = V (O)[B (O) + hC (O)], Sn = V(1)[B(1) + hC(1)], Tn С краевыми условиями умножим системы

T o —- V ( O ) B ( O ) ,           H o — V ( O ) f ( O ) ,

—-V(1)B(1), Hn — V(1)f (1).

(10 a ) , (10 b ) поступим следующим образом:

TT на матрицы D 0 , E 0

( D ol E o Г¹ I — G o ,

V ^x o J

и

(ENIFN ]f ' 1 — GN

V xN-1J и ENT, FNT соответственно. Получим:

D o ^T D o X 1 + D o ^T E o x o — D o ^T G o ,

E o ^T D o X 1 + E o T E o x o — E o T G o , TT  T

^E N ^E N x N ^{+ E} N F N x N - 1 ^E N ^G N , TT  T

F N ^E N x N ⁺ F N F N x N - 1 F N ^G N ^.

С краевыми условиями (11 a ) , (11 b ) поступим аналогично: умножим системы

( S o| T o

^{( S}n| T N

Ix'I I xo J и xN xN-1

= H o

= H N

на матрицы S ₀ ^T , T ₀ ^T и

S_N ^T , T_N ^T . Получим:

S o T S o X 1 + S o T Tx o — S o T H o ,

T O ^S O ^X 1 ⁺ T O T O x O ^— T O ^H o ,

S N S N x N ⁺ S N T N x N - 1 ^— S N H N ,

T N S N x N ⁺ T N T N x N - 1 ^— T N H N •

Итогом суммирования уравнений (12), (13), (16), (17) и (14), (15), (18), (19) является

K o ^x 1 + P O x o — q o ,

K N x N ⁺ P N x N - 1 ^{— q}N ,

( 2O )

где K _o, P _o, K N и P_N — квадратные ( n x n ) -матрицы, определенные по правилу

K о ^— D o D o + ^E 0 D o + S 0 S 0 + T o S o ,

P o ^— D o ^E o + ^E o ^E o + S o T O + T O T O ,

T   T   TT

K N ^N ^N ^{+ r} N ^N ⁺ S N S N ⁺ T N S N ,

TTTT

P N NN ^rN ^{+ r} N ^rN ⁺ S N T N ⁺ T N T N , q о = ^D о ^G о ⁺ E о ^G о ^{+ S} о H о ^{+ Т}о H о , q N = E N ^G N ⁺ F N ^G N ⁺ S N H N ⁺ T N H N . Для иллюстрации предлагаемого метода выбора недостающих краевых условий рассмотрим простой пример.

Пример 1. Система ко о

к

о I , t е [ о,1 ] . к ^e )

v                                  Ju (o) + u '(o) = 1,

Краевые условия имеют вид 1 _{z х z х}

( u ( 1 ) + u ' ( 1 ) = sin 1 + cos 1.

Точное решение данной системы u (t ) = sin t , v ( t ) = e* .

Применим описанный алгоритм поиска недостающих краевых усло- вий. Полагая A

: о

^о I = A + 0 )

(псевдообратная матрица), имеем

V (t)= 0 11

к 0 1)

. Умножим

исходную

систему слева на V ( t ) и вычислим

к

полученное выражение в точках t = 0 и t = 1.

В силу (5) получим:

: о о

Из второго уравнения видно, что

Обозначим x ( t ) = ( u ( t ), система (7) примет вид:

v ( t )У

отк e ¹ , t .

Jv (O)= 1, (v (1)=e.

и с учетом введенных ранее обозначений

к

г

(к

л

У 0 =

)

Л

IJ y1 =

к оJ

( sin 1 + cos 1 I ’

к

а системе (8) соответствует

к

(к

Л

, I ^y о =

(о I

к1)

( о I .

, J y1 =

к e

Далее применим двухточечную дискретизацию для исходных и полученных краевых условий и выполним дальнейшие указания описанного алгоритма. В результате, опустив все громоздкие выкладки, получим ис- комые краевые условия в виде:

к 0

(h + 1

к 0

л

x 1

+

(h — 1 0)

У

л

к

x N

—

У

1У

x 0 =

x N - 1

У

[ h'

к1У

( h (sin1 + cosOY к e У

Заключение

В статье не рассматривается конкретный вид трехточечной аппроксимации (9). Отметим, что ряд известных формул численного дифференци- рования нельзя применять для (1) даже с корректно заданными (согласо- ванными с правой частью) краевыми условиями

X ( 0 ) = x ₀

X (1) = Xn и в предположении, что решение задачи (1), (22) существует, единственно и достаточно гладкое. Это происходит, в частности, потому, что соотношения (9) и условия (22) дают вырожденную СЛАУ или процесс матричной прогонки является неустойчивым.

Также отметим, что предлагаемый подход не всегда дает СЛАУ (20), (21) полного ранга. Например, если применить предлагаемый подход для ДАУ (4), имеющей единственное решение, не зависящее от условий (3), и полагая a 1 (0) * 0, « 1 (1) * 0, то краевые условия (19), (20) имеют вид:

f X 1 (0) = ^ (0)

I X 1 (1) = ^ (1) .

Условия (3) для данного примера не заданы, а условий (23) недостаточно, чтобы формально реализовать метод матричной прогонки для данного примера (после трехточечной аппроксимации (9)).

В дальнейшем планируется выделить классы ДАУ (1) с краевыми условиями (3), для которых предлагаемый алгоритм выбора недостающих краевых условий и подходящая трехточечная аппроксимация порождали бы устойчивый процесс матричной прогонки.