Оценка порядка аппроксимации матричного метода численного интегрирования краевых задач для систем линейных неоднородных обыкновенных дифференциальных уравнений второго порядка с переменными коэффициентами. Сообщение 1. Краевые задачи с граничными условиями первого рода

Полный текст

Введение. Использование конечных разностей для аппроксимации производных в классическом методе сеток численного интегрирования краевых задач для неоднородных обыкновенных дифференциальных уравнений второго порядка (ОДУ2) с переменными коэффициентами с граничными условиями первого рода x + p(t)x + q(t)x = f (t), x(a) = x0 , x(b) = xn , (1) где p(t), q(t), f (t) - заданные функции, дифференцируемые нужное число раз, [a, b] - область интегрирования, x0 , xn - заданные числа, приводит ко второму порядку аппроксимации [1-6]. Такой же (второй) порядок аппроксимации имеют методы численного интегрирования ряда краевых задач для уравнений в частных производных [4-10]. Второй порядок обусловлен тем, что при аппроксимации производных конечными разностями было удержано всего три члена разложения в ряд Тейлора искомого решения задачи. Метод, использующий средства матричного исчисления, численного интегрирования разностных краевых задач для неоднородных ОДУ2, позволяющий увеличить количество удержанных членов до произвольного натурального числа в разложении искомой функции в ряд Тейлора, предложен в [11], при этом аппроксимация производных конечными разностями не использовалась. Оценка порядка аппроксимации предложенного в [11] матричного метода интегрирования краевых задач для неоднородных линейных ОДУ2 с граничными условиями первого, второго и третьего рода при различных значениях удержанного числа членов в рядах Тейлора дана в [12]. Исследуем возможность использования матричного метода для численного интегрирования систем неоднородных линейных ОДУ2 с переменными коэффициентами с различными граничными условиями на примере следующей системы для неизвестных функций x(t), y(t): u1 x + p1 x + r1 x + v1 y + q1 y + s1 y = f1 , u2 x + p2 x + r2 x + v2 y + q2 y + s2 y = f2 , (2) где uj , pj , rj , vj , qj , sj , fj - заданные функции аргумента t, дифференцируемые нужное число раз, j = 1, 2 - номер уравнения в системе (2), и поставим целью вычисление порядка аппроксимации. Далее будем придерживаться принятых в [4] обозначений: 1) D - область интегрирования, ограниченная отрезком [a, b], Dh - узлы сетки, определяемые значениями ti = t0 + ih, i = 1, 2, . . . , n, t0 = a, tn = b, h = (b - a)/n, n + 1 - число узлов сетки; 2) x(t), y(t) - непрерывные функции, являющиеся точным решением системы (2) с теми или иными граничными условиями; 390 Оценка порядка аппроксимации матричного метода численного интегрирования краевых . . . 3) [x]h , [y]h - сеточные функции, совпадающие с точным решением в узлах сетки Dh ; 4) x(h) , y (h) - искомые сеточные функции. Для краткости примем для любой функции обозначение ϕ(ti ) = ϕi , где ti - узел сетки Dh . В дальнейшем опустим индекс h в наименованиях сеточных функций [x]h , [y]h , x(h) , y (h) . 1. Краевые задачи с граничными условиями первого рода. Рассмотрим систему ОДУ2 (2) с граничными условиями первого рода x0 = x0 , y0 = y0 , xn = xn , yn = yn , (3) где x0 , y0 , xn , yn - заданные числа. В соответствии с матричным методом численного интегрирования выполним следующие преобразования. Для некоторых фиксированных k 2, k - степень используемого многочлена Тейлора, и номера узла сетки i, i = 1, 2, . . . , n - 1, запишем следующие многочлены Тейлора для функций x(t), y(t):  k h3 h2  k h x(k) ,  x - x + · · · + (-1) x = x - hx +  i-1 i i   2! i 3! i k! i   2 3  hk h h   yi-1 = yi - hyi + yi - yi + · · · + (-1)k yi(k) , 2! 3! k! (4)  h2 h3 hk (k)   xi+1 = xi + hxi + xi + xi + · · · + xi ,   2! 3! k!    2 3 k  h h h  yi+1 = yi + hy + y + y + · · · + y (k) , i i i 2! 3! k! i вычислим производные по аргументу t от обеих частей уравнений системы (2) и запишем их в узле ti в виде  (u1i x + p1i x + r1i xi + v1i y + q1i y + s1i yi )(r) = (f1i )(r) , i i i i (5) (u x + p x + r x + v y + q y + s y )(r) = (f )(r) , 2i i 2i i 2i i 2i i 2i i 2i i 2i где r = 1, 2, . . . , k - 2; здесь и далее в наименованиях функций первый из пары нижних индексов означает номер уравнения в системе (2), второй - номер узла сетки. Из многочленов Тейлора (4), дифференциальных уравнений системы (2), записанных в узле ti , и производных (5) составим следующую систему из 2k + 2 уравнений:  k h2 h3  k h x(k) = x  x - hx + x - x + · · · + (-1)  i i-1 , i   2! i 3! i k! i   k  h3 h2  k h y (k) = y   y - y + · · · + (-1) y - hy + i i-1 , i   2! i 3! i k! i    h2 h3 hk (k)    xi + hxi + xi + xi + · · · + x = xi+1 , (6)   2! 3! k! i    h2 h3 hk (k)    y + y + · · · + y = yi+1 , y + hy + i  i i i  2! 3! k! i    r1i xi + p1i xi + u1i xi + s1i yi + q1i yi + v1i yi = f1i ,     r2i xi + p2i xi + u2i xi + s2i yi + q2i yi + v2i yi = f2i , 391 М а к л а к о в В. Н.  r1i xi + (r1i + p1i )xi + (p1i + u1i )xi + u1i xi +    +s1i yi + (s1i + q1i )yi + (q1i + v1i )yi + v1i yi = f1i ,       r2i xi + (r2i + p2i )xi + (p2i + u2i )xi + u2i xi +    +s2i yi + (s2i + q2i )yi + (q2i + v2i )yi + v2i yi = f2i ,      ..    .    (k-2) (k) (k-2) (k) (k-2)   r1i xi + · · · + u1i xi + s1i yi + · · · + v1i yi = f1i ,     (k-2) (k) (k-2) (k) (k-2) r2i xi + · · · + u2i xi + s2i yi + · · · + v2i yi = f2i . Система (6) - система линейных алгебраических уравнений (СЛАУ), имеющая в матричной форме вид Aki W ki = Gki , где  Aki 1 0 1 0 r1i r2i ... 0 1 0 1 s1i s2i ...        =      (k-2) (k-2) r1i s1i (k-2) (k-2) r2i s2i 2 h 0 -h 0 2! h2 0 -h 0 2! h2 h 0 0 2! 2 0 h 0 h2! p1i q1i u1i v1i p2i q2i u2i v2i ... ... ... ... ... ... ... ... ... ... ... ... k . . . (-1)k hk! ... 0 hk ... k! ... 0 ... 0 ... 0 ... ... ... u1i ... u2i (k) W ki = xi yi xi yi xi yi . . . xi Gki = xi-1 yi-1 xi+1 yi+1 f1i f2i . . . (k) yi (k-2) f1i  0 k (-1)k hk!     0  hk  k!  , 0   0  ...   v1i  v2i (7) , (k-2) f2i . Здесь и далее верхний индекс k означает степень используемого многочлена Тейлора, если речь не идет о показателях алгебраических степеней и степенях производных; второй из пары верхних индексов (если он существует) в наименованиях матриц и их элементов означает номер узла сетки. Матрицу Aki будем называть локальной матрицей. Предполагая существование обратной матрицы B ki = (Aki )-1 от локальной матрицы (7), запишем матричное равенство B ki Gki = W ki , или следующие соотношения в координатной форме: ki ki ki ki ki bki 11 xi-1 + b12 yi-1 + b13 xi+1 + b14 yi+1 + b15 f1i + b16 f2i + · · · + (k-2) + bki 1 (2k+1) f1i (k-2) = xi , (8) (k-2) = yi , (9) + bki 1 (2k+2) f2i ki ki ki ki ki bki 21 xi-1 + b22 yi-1 + b23 xi+1 + b24 yi+1 + b25 f1i + b26 f2i + · · · + (k-2) + bki 2 (2k+1) f1i 392 + bki 2 (2k+2) f2i Оценка порядка аппроксимации матричного метода численного интегрирования краевых . . . ki ki ki ki ki bki 31 xi-1 + b32 yi-1 + b33 xi+1 + b34 yi+1 + b35 f1i + b36 f2i + · · · + (k-2) + bki 3 (2k+1) f1i (k-2) = xi , (10) (k-2) = yi , (11) + bki 3 (2k+2) f2i ki ki ki ki ki bki 41 xi-1 + b42 yi-1 + b43 xi+1 + b44 yi+1 + b45 f1i + b46 f2i + · · · + (k-2) + bki 4 (2k+1) f1i + bki 4 (2k+2) f2i ... ki ki ki ki bki (2k+1) 1 xi-1 + b(2k+1) 2 yi-1 + b(2k+1) 3 xi+1 + b(2k+1) 4 yi+1 + b(2k+1) 5 f1i + (k-2) ki + bki (2k+1) 6 f2i + · · · + b(2k+1) (2k+1) f1i (k-2) + bki (2k+1) (2k+2) f2i (k) = xi , (12) ki ki ki ki bki (2k+2) 1 xi-1 + b(2k+2) 2 yi-1 + b(2k+2) 3 xi+1 + b(2k+2) 4 yi+1 + b(2k+2) 5 f1i + (k-2) ki + bki (2k+2) 6 f2i + · · · + b(2k+2) (2k+1) f1i (k-2) + bki (2k+2) (2k+2) f2i (k) = yi , (13) ki где bki jm - элементы матрицы B . Из соотношений (8), (9), являющихся разностными уравнениями второго порядка для трехточечного шаблона ti-1 , ti , ti+1 , i = 1, 2, . . . , n - 1, и связывающих значения xi-1 , yi-1 , xi , yi , xi+1 , yi+1 искомых функций, составим, учитывая граничные условия (3), систему  k1 x - bk1 y = bk1 x + bk1 y + x1 - b13  2 14 2 11 0 12 0   (k-2) (k-2)  k1 f + bk1 f + · · · + bk1  +b + bk1 ,  15 11 16 21 1 (2k+1) f11 1 (2k+2) f21    k1 k1 k1   y1 - bk1 23 x2 - b24 y2 = b21 x0 + b22 y0 +    (k-2) (k-2) k1 k1  +bk1 + bk1 ,  25 f11 + b26 f21 + · · · + b2 (2k+1) f11  2 (2k+2) f21   .   ..     ki ki ki  -bki  11 xi-1 - b12 yi-1 + xi - b13 xi+1 - b14 yi+1 =   (k-2) (k-2)  ki ki  = bki + bki , 15 f1i + b16 f2i + · · · + b1 (2k+1) f1i 1 (2k+2) f2i (14) ki ki ki ki  -b21 xi-1 - b22 yi-1 + yi - b23 xi+1 - b24 yi+1 =    (k-2) (k-2) ki ki  = bki + bki ,  25 f1i + b26 f2i + · · · + b2 (2k+1) f1i  2 (2k+2) f2i   .   ..     k (n-1) k (n-1) k (n-1) k (n-1)   -b11 xn-2 - b12 yn-2 + xn-1 = b13 xn + b14 yn +    k (n-1) k (n-1) k (n-1) (k-2) k (n-1) (k-2)   +b15 f1 (n-1) + b16 f2 (n-1) + · · · + b1 (2k+1) f1 (n-1) + b1 (2k+2) f2 (n-1) ,    k (n-1) k (n-1) k (n-1) k (n-1)   -b21 xn-2 - b22 yn-2 + yn-1 = b23 xn + b24 yn +    k (n-1) k (n-1) (k-2) k (n-1) (k-2)  +bk (n-1) f +b f + ··· + b f +b f . 25 1 (n-1) 26 2 (n-1) 2 (2k+1) 1 (n-1) 2 (2k+2) 2 (n-1) 393 М а к л а к о в В. Н. Отметим, что в центральном узле шаблона уравнения системы (14) с нечетными номерами содержат искомые значения xi , с четными номерами - значения yi . Система (14) является СЛАУ, имеющей блочную трехдиагональную матрицу, решение которой позволяет найти метод матричной прогонки [5, 13], что дает возможность вычислить искомые значения xi , yi во внутренних узлах (r) (r) сетки ti , i = 1, 2, . . . , n - 1. Производные xi , yi , r = 1, 2, . . . , k, во внутренних узлах сетки найдем по формулам (10)-(13) с использованием найденных значений xi , yi . Фактически СЛАУ (14) есть запись в развернутом виде разностной краевой задачи для системы двух разностных уравнений второго порядка с граничными условиями первого рода для отыскания вектор-функции zi = [xi yi ] , i = 1, 2, . . . , n - 1, удовлетворяющей следующим равенствам [4, 14]: 1 Bjki zi+j = g ki , i = 1, 2, . . . , n - 1, λ0 z0 = z0 , λn zn = zn , (15) j=-1 где ki B-1 = g ki = ki -bki 11 -b12 , ki -b21 -bki 22 g1i = C ki F ki , g2i 1 0 , 0 1 B0ki = λ0 = λn = B1ki = ki -bki 13 -b14 , ki -b23 -bki 24 ki ki ki ki ki bki 15 b16 b17 b18 . . . b1 (2k+1) b1 (2k+2) ki ki ki i ki bki 25 b26 b27 b28 . . . b2 (2k+1) b2 (2k+2) C ki = (k-2) F ki = f1i f2i f1i f2i . . . f1i (k-2) f2i , z0 = [x0 y0 ] , , zn = [xn yn ] . 2. Предварительные оценки. При некоторых фиксированных номере k 2 и номере i, i = 1, 2, . . . , n - 1, запишем следующие ряды Тейлора для функции x(t): xi-1 = xi - hxi + h3 h2 xi - x + ··· = 2! 3! i k = mh (-1) m=0 xi+1 = xi + hxi + m m! ∞ (m) xi (-1)m + m=k+1 hm (m) k k x = Px,i-1 + Rx,i-1 , m! i h2 h3 xi + x + ··· = 2! 3! i k = m=0 hm (m) x + m! i ∞ m=k+1 hm (m) k k x = Px,i+1 + Rx,i+1 , m! i где k k Rx,i-1 , Rx,i+1 = 394 hk+1 (k+1) x (ξi ) = O(hk+1 ), (k + 1)! ξi ∈ (ti , ti+1 ) Оценка порядка аппроксимации матричного метода численного интегрирования краевых . . . - дополнительные члены разложений в ряд Тейлора в форме Лагранжа [15]. Из двух последних точных равенств найдем ∞ k k Rx,i+1 - Rx,i-1 = m=k+1 hm (m) x - m! i ∞ (-1)m m=k+1 hm (m) x = m! i ∞ = m=k+1 hm (m) x (1 - (-1)m ) (16) m! i и, аналогично, ∞ k k Rx,i+1 + Rx,i-1 = m=k+1 hm (m) x + m! i ∞ (-1)m m=k+1 hm (m) x = m! i ∞ = m=k+1 hm (m) x (1 + (-1)m ) . (17) m! i Из равенств (16), (17) для четного k имеем k k Rx,i+1 - Rx,i-1 = k k Rx,i+1 + Rx,i-1 = hk+1 (k+1) x (1 + 1) + (k + 1)! i hk+2 (k+2) x (1 - 1) + · · · = O(hk+1 ), (18) + (k + 2)! i hk+1 (k+1) x (1 - 1) + (k + 1)! i hk+2 (k+2) + xi (1 + 1) + · · · = O(hk+2 ) (19) (k + 2)! и для нечетного k - k k Rx,i+1 - Rx,i-1 = k k Rx,i+1 + Rx,i-1 = hk+1 (k+1) x (1 - 1) + (k + 1)! i hk+2 (k+2) + xi (1 + 1) + · · · = O(hk+2 ), (20) (k + 2)! hk+1 (k+1) x (1 + 1) + (k + 1)! i hk+2 (k+2) + x (1 - 1) + · · · = O(hk+1 ). (21) (k + 2)! i Замена в равенствах (18)-(21) x на y даст аналогичные оценки для доk k полнительных членов Ry,i-1 , Ry,i+1 . 395 М а к л а к о в В. Н. 3. Оценка порядка аппроксимации разностной краевой задачи с граничными условиями первого рода. Соответствующую дифференциальной краевой задаче (1) разностную краевую задачу будем обозначать, по аналогии с [4], в компактной символической форме как Lkh x = fhk . При исследовании дифференциальной краевой задачи (1) для одного уравнения сеточная функция xi , i = 0, 1, . . . , n, являющаяся решением разностной краевой задачи, при подстановке в уравнения этой разностной краевой задачи обратит их в верные равенства. В [4] показано, что подстановка в уравнения задачи сеточной функции [xi ], отличающейся от xi , приведет к некоторому отличию от верных равенств. Эти отличия и характеризует невязка δfhk [4]. Иными словами, подстановка [x] в Lkh x = fhk приведет к Lkh [x] = fhk + δfhk . В работе [4] в качестве оценки величины невязки принята норма δfhk = max |δfhk0 |, |δfhkn |, |δfhki | , i = 1, 2, . . . , n - 1, где первые две компоненты характеризуют меру отличий в граничных узлах сетки Dh , оставшиеся - во внутренних узлах. Согласно [4, 6] разностная краевая задача (14) аппроксимирует дифференциальную краевую задачу (1) на точном решении x, если δfhk → 0 при h → 0. Если при этом имеет место неравенство δfhk Chk , где C > 0, k > 0 - некоторые постоянные, не зависящие от h, то говорят, что имеет место аппроксимация порядка k относительно величины h. Вычислим порядок аппроксимации рассматриваемой краевой задачи. СЛАУ (14), или, что то же самое, разностную краевую задачу (15) запишем, по аналогии с [4], в компактной символической форме как Lkh z = fhk , или, в дальнейшем для краткости, как Lkh , где Lkh z ≡ Lkh 396 x y  ki b   xi-1 - - 11   bki  15     bki  21  - ki xi-1 - b 26 =  x0 ,      y0 ,     xn ,    yn , bki xi bki bki 12 13 14 y + - x - yi+1 , i-1 i+1 ki ki ki b15 b15 b15 bki 15 bki yi bki bki 22 23 24 y + - x - y , i-1 i+1 ki ki ki i+1 bki b b b 26 26 26 26 (22) Оценка порядка аппроксимации матричного метода численного интегрирования краевых . . .   bki  16  f + f +  1i  ki 2i  b  15      bki   25 f1i + f2i + ki fhk = b26   x0 ,     y0 ,      xn ,    yn , bki bki bki bki 1 (2k+1) (k-2) 1 (2k+2) (k-2) 17 18 + , f f f f2i + + · · · + 1i 2i 1i ki ki ki b15 b15 b15 bki 15 bki bki bki bki 2 (2k+2) (k-2) 2 (2k+1) (k-2) 27 28 f1i + ki f2i + · · · + f1i + f2i , ki ki b26 b26 b26 bki 26 для всех i = 1, 2, . . . , n - 1. ki Как будет показано ниже, выбор того или иного знаменателя (bki 15 или b26 ) в двух первых равенствах задачи (22) не является существенным при оценке порядка аппроксимации. Для задачи (22) обобщение приведенного выше определения порядка аппроксимации запишем в форме Lkh [z] = fhk + δfhk , (23) где k k δfhk = δf1h δf2h ; k , δf k - величины невязок в узлах сетки D , появление которых обусловδf1h h 2h лено соответственно первым и вторым уравнениями задачи (22). Итоговую норму невязку задачи в соответствии с [4] определим как k k δfhk = max δf1h , δf2h . (24) При фиксированных k, i, i = 1, 2, . . . , n - 1, имеем точные равенства  k h2 h3  k h [x(k) ] = [x k  [x ] - h[x ] + [x ] - [x ] + · · · + (-1)  i i-1 ] - Rx,i-1 , i   2! i 3! i k! i    h2 h3 hk   k [yi ] - h[yi ] + [yi ] - [yi ] + · · · + (-1)k [yi(k) ] = [yi-1 ] - Ry,i-1 , 2! 3! k! (25)  h2 h3 hk (k)  k  [xi ] + h[xi ] + [xi ] + [xi ] + · · · + [xi ] = [xi+1 ] - Rx,i+1 ,   2! 3! k!    2 3 k  h h h  (k) k  [yi ] + h[yi ] + [yi ] + [yi ] + · · · + [yi ] = [yi+1 ] - Ry,i+1 . 2! 3! k! Составим СЛАУ, в которую внесем равенства (25), дифференциальные уравнения системы (2), записанные в узле ti , и точные равенства (5). В итоге получим 397 М а к л а к о в В. Н.  k h3 h2  k k h [x(k) ] = [x  [x [x ] - ] + · · · + (-1) [x ] - h[x ] +  i-1 ] - Rx,i-1 , i i i i   2! 3! k! i    h3 h2 hk (k)   k  , [yi ] - h[yi ] + [yi ] - [yi ] + · · · + (-1)k [yi ] = [yi-1 ] - Ry,i-1   2! 3! k!    2 3 k  h h h (k)  k  [x ] + h[x ] + [x ] + [x ] + · · · + [xi ] = [xi+1 ] - Rx,i+1 ,  i i i i  2! 3! k!     h3 hk (k) h2  k   [y [y [yi ] = [yi+1 ] - Ry,i+1 [y ] + h[y ] + ] + · · · + , ] + i  i i i  2! 3! k!     r1i [xi ] + p1i [xi ] + u1i [xi ] + s1i [yi ] + q1i [yi ] + v1i [y ] = f1i ,    r2i [xi ] + p2i [xi ] + u2i [xi ] + s2i [yi ] + q2i [yi ] + v2i [yi ] = f2i ,     r1i [xi ] + (r1i + p1i )[xi ] + (p1i + u1i )[xi ] + u1i [xi ]+      +s1i [yi ] + (s1i + q1i )[yi ] + (q1i + v1i )[yi ] + v1i [yi ] = f1i ,       r2i [xi ] + (r2i + p2i )[xi ] + (p2i + u2i )[xi ] + u2i [xi ]+     +s2i [yi ] + (s2i + q2i )[yi ] + (q2i + v2i )[yi ] + v2i [yi ] = f2i ,       ..   .    (k-2) (k) (k-2) (k) (k-2)   r1i [xi ] + · · · + u1i [xi ] + s1i [yi ] + · · · + v1i [yi ] = f1i ,      (k-2) (k) (k-2) (k) (k-2) r2i [xi ] + · · · + u2i [xi ] + s2i [yi ] + · · · + v2i [yi ] = f2i . (26) В матричной форме система (26) имеет вид Aki [W ki ] = [Gki ] , где  [xi ]  [yi ]   [x ]   i   [y ]   i   [x ]   i    [W ki ] =  [yi ]  ,  [xi ]     [yi ]     ...   (k)  [xi ]  (k) [yi ]   k [xi-1 ] - Rx,i-1 k   [yi-1 ] - Ry,i-1    [xi+1 ] - Rk x,i+1     [y ] - Rk  i+1 y,i+1    f1i   ki , [G ] =  f2i     f 1i     f 2i     ...   (k-2)   f 1i (k-2) f2i а матрица Aki определяется формулой (7). В предположении существования обратной матрицы B ki = (Aki )-1 от матрицы Aki найдем B ki [Gki ] = [W ki ]. Выпишем первые два уравнения последнего матричного равенства: 398 Оценка порядка аппроксимации матричного метода численного интегрирования краевых . . . k ki k ki k bki 11 ([xi-1 ] - Rx,i-1 ) + b12 ([yi-1 ] - Ry,i-1 ) + b13 ([xi+1 ] - Rx,i+1 )+ k ki ki ki ki + bki 14 ([yi+1 ] - Ry,i+1 ) + b15 f1i + b16 f2i + b17 f1i + b18 f2i + (k-2) + · · · + bki 1 (2k+1) f1i (k-2) + bki 1 (2k+2) f2i = [xi ], k ki k ki k bki 21 ([xi-1 ] - Rx,i-1 ) + b22 ([yi-1 ] - Ry,i-1 ) + b23 ([xi+1 ] - Rx,i+1 )+ k ki ki ki ki + bki 24 ([yi+1 ] - Ry,i+1 ) + b25 f1i + b26 f2i + b27 f1i + b28 f2i + (k-2) + · · · + bki 2 (2k+1) f1i (k-2) + bki 2 (2k+2) f2i = [yi ], ki в узле t , или, после преобразований, где bki i jm - элементы матрицы B - bki [xi ] bki bki bki 12 13 14 11 [x ] - [y ] + - [x ] - [yi+1 ] = f1i + i-1 i-1 i+1 ki ki ki ki b15 b15 b15 b15 bki 15 + bki bki bki bki bki 1 (2k+1) (k-2) 1 (2k+2) (k-2) 16 17 18 f + f + f + · · · + f + f2i - 2i 1i ki 1i ki 2i ki bki b b b bki 15 15 15 15 15 - - k ki k bki 11 Rx,i-1 + b13 Rx,i+1 bki 15 - k ki k bki 12 Ry,i-1 + b14 Ry,i+1 bki 15 , (27) bki [yi ] bki bki bki bki 22 23 24 25 21 [x ] - [y ] + - [x ] - [y ] = f1i + i-1 i-1 i+1 i+1 ki ki ki ki ki ki b26 b26 b26 b26 b26 b26 + f2i + bki bki bki bki 2 (2k+1) (k-2) 2 (2k+2) (k-2) 27 28 f + f + · · · + f + f2i - 1i 1i ki 2i ki bki b b bki 26 26 26 26 - k ki k bki 21 Rx,i-1 + b23 Rx,i+1 bki 26 - k ki k bki 22 Ry,i-1 + b24 Ry,i+1 bki 26 . (28) Отбрасывание двух последних дробей в составленной из уравнений (27), (28) во внутренних узлах сетки Dh системе, что равносильно переходу от точного решения [xi ], [yi ] к искомому приближенному xi , yi , приведет эту систему к разностной краевой задаче (22). Следовательно, в соответствии с (23) две последние дроби в уравнениях (27), (28) характеризуют величину невязки в узлах ti , i = 1, 2, . . . , n - 1; в итоге для рассматриваемой задачи имеем k ki k k ki k bki bki 11 Rx,i-1 + b13 Rx,i+1 12 Ry,i-1 + b14 Ry,i+1 ki δf1h =- - , (29) bki bki 15 15 ki δf2h =- k ki k bki 21 Rx,i-1 + b23 Rx,i+1 bki 26 - k ki k bki 22 Ry,i-1 + b24 Ry,i+1 bki 26 . (30) Первое слагаемое в равенствах (29), (30) характеризует величину невязки, появление которой обусловлено функцией x, второе - функцией y. Для краевой задачи с граничными условиями первого рода величины невязок в граничных узлах сетки обращаются в нуль в силу того, что два 399 М а к л а к о в В. Н. первых уравнения задачи (22) при i = 1, n - 1 содержат заданные значения x0 , y0 , xn , yn искомых функций. ki и δf ki во внутренних узлах сетки D . Исследуем невязки δf1h h 2h В узле ti введем компактные обозначения для определителей второго порядка, например u v (31) Ui Vi = 1i 1i , u2i v2i где, очевидно, Ui Vi = -Vi Ui . В дальнейшем будем опускать индекс i в обозначениях определителей вида (31). Непосредственными вычислениями можно для любого k 3 убедиться в справедливости оценки ki 2i M11 ≈ (U V )k-2 M11 , (32) ki - алгебраическое дополнение элемента aki транспонированной логде M11 11 кальной матрицы Aki . Отсутствие в настоящей работе полного вывода формулы (32) обусловлено лишь громоздкостью выкладок. Для краевых задач для одного ОДУ2 полный вывод формулы (32) дан в [12]. Тем не менее рассмотрим частный случай при k = 3, опуская индекс i - номер узла сетки. Воспользуемся в дальнейшем известными свойствами определителя [16]. Имеем, пренебрегая старшими степенями, следующее: 3 M11 1 0 -h = 02 h 2 0 h 0 1 0 h 0 h2 2 0 h3 3! 0 3 - h3! 0 1 0 0 h -h 0 + u1 0 h2 2 2 h 0 23 - h3! 0 s1 s2 s1 s2 p1 p2 p1 p2 q1 q2 q1 q2 u1 u2 u1 u2 = h2 v1 v2 v1 v2 2 0 0 0 u1 u2 h3 0 0 v1 v2 3! 1 s1 s2 s2 0 p1 p2 p2 h q1 q2 q2 0 u1 u2 u2 - u2 h2 v1 v2 v2 2 h3 0 0 v2 3! 1 0 3 -h h 3! 02 h 23 - h3! 1 0 -h 0 h2 23 - h3! 0 h 0 h2 2 0 0 1 0 h 0 h2 2 h3 3! 1 0 h 0 h2 2 h3 3! s1 s2 s1 s2 p1 p2 p1 p2 q1 q2 q1 q2 u1 u2 u1 u2 + v1 v2 v1 v2 0 0 v1 v2 s1 s2 s1 p1 p2 p1 q1 q2 q1 u1 u2 u1 = v1 v2 v1 0 0 v1 2 = c3 h6 + O h7 + u1 d3 h4 + O h5 + v2 M11 - 2 2 2 - u2 g3 h4 + O h5 + v1 M11 ≈ (u1 v2 - u2 v1 ) M11 = U V M11 , где pj , qj , uj , vj , j = 1, 2, - функции от pj , qj , uj , vj и их первых производных, c3 , d3 , g3 - не зависящие от h величины. Формулы, аналогичные (32), имеют место, по крайней мере, для первых шести элементов первой и второй строк матрицы (Aki ) ; на основании чего и очевидных равенств bki 1j bki 15 400 = k M1j k M15 , bki 2j bki 26 = k M2j k M26 , j = 1, 2, 3, 4, Оценка порядка аппроксимации матричного метода численного интегрирования краевых . . . следуют оценки bki 1j bki 15 ≈ 2 M1j 2 M15 , bki 2j bki 26 ≈ 2 M2j 2 M26 , j = 1, 2, 3, 4. (33) Невязка (29) с учетом соотношений (33) примет вид ki δf1h ≈- 2 Rk 2 k M11 x,i-1 + M13 Rx,i+1 2 M15 - 2 Rk 2 k M12 y,i-1 + M14 Ry,i+1 2 M15 . (34) Вычислим оценки первых пяти алгебраических дополнений первой строки матрицы A2i , i = 1, 2, . . . , n - 1. Имеем, пренебрегая старшими степенями, следующее: 2 M11 1 0 = -h 0 0 h 0 1 0 h 0 h2 2 h2 h2 2 0 2 s1 p1 q1 u1 v1  s2 0 p2 1 q2 = -2 h2 h 2 u2 0 v2 1 2 = -2 h - 0 h2 2 1 0 0 h2 2 s1 s2 p1 p2 u1 u2 = v1 v2 1 s1 s2 h 0 u1 u2 + 2 h2 v1 v2 2  s1 s2 p1 p2  = v1 v2 h h2 h2 SU + PV + SP = 2 2 2 h h2 h3 U V - P V + SU - SP ≈ h2 (2 U V - h P V ) , (35) 2 2 4 = -2 h2 -U V - = 2 h2 где были использованы обозначения определителей, аналогичные введенным формулой (31), 2 M13 0 1 -h 0 = 02 -h h 0 2 h2 0 2 1 0 h 0 s1 s2 p1 p2 q 1 q2 = u1 u2 v1 v2 h2 2 = 2 h2 U V + 2 M12 0 -h = - 02 0 h 0 h 2 h2 2 0 0 1 0 h 0 h2 2 h h2 h3 P V + SU + SP 2 2 4 s1 s2 p1 p2 q1 q2 = -h3 u1 u2 v1 v2 0 0 1 0 1 h 0 h2 2 ≈ h2 (2 U V + h P V ) , (36) s1 s2 q1 q2 u1 u2 = v1 v2 401 М а к л а к о в В. Н. = -h 3 1 h h2 2 s1 s2 h2 q1 q2 = h3 -QV + h SV - SQ ≈ h3 (-QV + h SV ) , (37) 2 v1 v2 2 M14 2 M15 0 1 -h 0 = 02 -h h 0 2 h2 0 2 0 1 -h 0 = - 02 -h h 0 2 h2 0 2 0 h 0 h2 2 0 1 0 h 0 h2 2 0 h 0 h2 2 0 s1 s2 p 1 p2 q1 q2 ≈ h3 (QV + h SV ) , u1 u2 v1 v2 s2 p2 q2 = -4 h3 u2 v2 0 = 2 h4 1 0 1 0 h2 2 0 0 1 0 0 h2 0 0 1 h 0 h2 2 (38) s2 q2 u2 = v2 s2 h2 u2 = 2 h4 -v2 + s2 2 v2 ≈ -2 h4 v2 . (39) С учетом оценок (35)-(39) величину невязки (34), появление которой обусловлено первым уравнением задачи (22), на точном решении [xi ], [yi ] во внутренних узлах ti сетки Dh задачи Lkh запишем как ki δf1h ≈- k k h2 (2 U V - h P V ) Rx,i-1 + (2 U V + h P V ) Rx,i+1 - -2 h4 v2 k k h3 (-QV + h SV ) Ry,i-1 + (QV + h SV ) Ry,i+1 - = -2 h4 v2 k k k k 2 U V (Rx,i+1 + Rx,i-1 ) + h P V (Rx,i+1 - Rx,i-1 ) = + 2 h2 v2 k k k k QV (Ry,i+1 - Ry,i-1 ) + h SV (Ry,i+1 + Ry,i-1 ) + 2 h v2 (40) для произвольного k 2. Из оценки (40) для четного k с учетом (18), (19) имеем ki δf1h ≈ 2 U V O(hk+2 ) + h P V O(hk+1 ) QV O(hk+1 ) + h SV O(hk+2 ) + = 2 h2 v2 2 h v2 = O(hk ) + O(hk ) + O(hk ) + O(hk+2 ) = O(hk ), (41) а для нечетного k с учетом (20), (21) имеем ki δf1h ≈ 402 2 U V O(hk+1 ) + h P V O(hk+2 ) QV O(hk+2 ) + h SV O(hk+1 ) + = 2 h2 v2 2 h v2 Оценка порядка аппроксимации матричного метода численного интегрирования краевых . . . = O(hk-1 ) + O(hk+1 ) + O(hk+1 ) + O(hk+1 ) = O(hk-1 ) (42) для всех i = 1, 2, . . . , n - 1. Исследование невязки (30) дало схожий результат: ki δf2h ≈ O(hk ) (43) ki δf2h ≈ O(hk-1 ) (44) для четного k и для нечетного k. Из оценки (41) следует, что при четном k вклады функций x, y в величиki сравнимы между собой, тогда как оценка (42) показывает, ну невязки δf1h что при нечетном k обусловленная функцией x составляющая невязки как имеющая более низкий порядок играет главенствующую роль по сравнению с составляющей, обусловленной функцией y. Практически аналогичной окаki с тем лишь отличием, что при нечетном k залась ситуация с невязкой δf2h главенствующую роль во вкладе в величину невязки играет функция y. В итоге в соответствии с (24) норму невязки в равенстве (23) запишем для четного k с учетом (41), (43) как δfhk = max k δf1h k δf2h = max k0 |, |δf k1 |, . . . , |δf kn | max |δf1h 1h 1h k0 k1 |, . . . , |δf kn | max |δf2h |, |δf2h 2h = max O(hk ) O(hk ) = = O(hk ), (45) и для нечетного k с учетом (42), (44) как δfhk = O(hk-1 ). (46) 2m+1 , где Анализ норм невязок (45), (46) показывает, что задачи L2m h , Lh m - натуральное число, имеют одинаковый порядок аппроксимации. Следовательно, на практике для уменьшения объема вычислений следует использовать задачи Lkh с четными значениями k. Действительно, число арифметических операций только для нахождения обратной матрицы от матрицы 2 размерности n методом Гаусса вычисляется по формуле 38 n3 - n2 - n6 - 1 [17]. Вывод о правомерности выбора четного значения k ниже будет подтвержден численным экспериментом. 2 ≈ -2h4 u имеет, как следует Значение алгебраического дополнения M26 1 2 , из (39), одинаковый порядок малости с алгебраическим дополнением M15 2 2 что свидетельствует о несущественной роли в использовании M15 или M26 ki или, что то же самое, в использовании bki 15 или b26 при вычислении невязок по формулам (29), (30). 4. Оценка погрешностей. При выполнении численного эксперимента в качестве суммарной оценки относительных погрешностей использованы следующие нормы: Dxk = n 2 i=0 (xi - [xi ]) n i=0 [xi ] · 100 %, Dyk = n 2 i=0 (yi - [yi ]) n i=0 [yi ] · 100 %, (47) 403 М а к л а к о в В. Н. которые можно трактовать как некий аналог коэффициента вариации в статистике, характеризующий меру разброса в процентах [18], а в качестве максимальной оценки абсолютных погрешностей: Exk = max xi - [xi ] , Eyk = max yi - [yi ] , i = 0, 1, . . . , n. (48) В качестве примера использована имеющая аналитическое решение система нелинейных ОДУ2 вида (2)  t2 + 2   x - tx - x - ty = 2t cos t, t2 (49)  1 x t2 - 4   x + +y + y = -2t sin t 2 t 2t2 с граничными условиями x(2π) = 0, x(3π) = 0, y(2π) = 4π 2 , y(3π) = -9π 2 . (50) В вычислениях использовались следующие параметры сетки: n = 15, h = 0.20944. Результаты численного эксперимента для краевой задачи (49) с граничными условиями (50) приведены в табл. 1, 2. В табл. 2 нормы Dxk , Dyk , Exk , Eyk для производных x (t), y (t) характеризуют суммарные оценки относительных погрешностей. Максимальные оценки абсолютных погрешностей, соответственно, вычислены по формулам (47), (48), в которых значения функций заменены на значения своих первых производных, найденных по формулам (10), (11), во внутренних узлах обла(j) (j) сти интегрирования; найденных при вычисленных x1 , y1 , j = 1, 2, . . . , k, по (10)-(13) при i = 1 по формулам h2 hk-1 (k) x1 - · · · + (-1)k-1 x , 2! (k - 1)! 1 h2 hk-1 (k) y0 = y1 - hy1 + y1 - · · · + (-1)k-1 y 2! (k - 1)! 1 x0 = x1 - hx1 + (51) (52) в левой границе t0 = a; найденных xn , yn по аналогичным (51), (52) разложениям в правой границе tn = b. Таблица 1 Значения погрешностей для решения краевой задачи (49), (50) [The values of the errors for the solution of the boundary value problem (49), (50)] k 2 3 4 5 6 7 Dxk , % Dyk , % Exk Eyk 404 1.62 · 10-1 1.21 · 10-1 3.59 · 10-1 2.19 · 10-1 1.47 · 10-1 2.13 · 10-1 3.47 · 10-1 3.56 · 10-1 4.21 · 10-4 6.96 · 10-4 1.01 · 10-3 1.15 · 10-3 1.70 · 10-4 4.49 · 10-4 4.34 · 10-4 7.21 · 10-4 1.57 · 10-5 9.85 · 10-6 3.45 · 10-5 1.77 · 10-5 1.60 · 10-5 1.09 · 10-5 3.52 · 10-5 1.97 · 10-5 Оценка порядка аппроксимации матричного метода численного интегрирования краевых . . . Таблица 2 Значения погрешностей для первых производных решения краевой задачи (49), (50) [The values of the errors for the first derivatives of the boundary value problem (49), (50)] k 2 3 4 5 6 7 Dxk , % Dyk , % Exk Eyk 4.53 · 10-1 3.43 · 10-1 7.32 · 10-1 8.20 · 10-1 1.98 · 10-1 1.50 · 10-1 4.07 · 10-1 3.90 · 10-1 1.31 · 10-3 1.09 · 10-3 2.27 · 10-3 2.57 · 10-3 2.80 · 10-4 2.90 · 10-4 5.09 · 10-4 6.94 · 10-4 3.78 · 10-5 2.34 · 10-5 7.09 · 10-5 5.66 · 10-5 3.86 · 10-5 2.41 · 10-5 7.08 · 10-5 5.88 · 10-5 Практически аналогичный характер изменения погрешностей (динамика и абсолютные значения) имел место для ряда систем ОДУ2, в частности для системы (1 + t)x + 2x + ty - 2y = 2 sin 2t, x + 2x - 2ty = 2(1 + t2 ) sin 2t с граничными условиями x(2π) = 0, x(3π) = 0, y(2π) = 2π, y(3π) = 3π. Анализ данных табл. 1, 2 свидетельствует, что суммарные относительные Dxk , Dyk , Dxk , Dyk и максимальные абсолютные Exk , Eyk , Exk , Eyk погрешности 2m+1 , имеющих одинаковый порядок аппроксимации, различазадач L2m h и Lh ются незначительно для любого натурального m = 1, 2, 3. Выводы. 1. Теоретически выявлены закономерности между порядком аппроксимации матричного метода и степенью k используемого многочлена Тейлора в разностных краевых задачах для систем линейных ОДУ2 с граничными условиями первого рода. Установлено, что порядок аппроксима2m+1 совпадает и равен 2m для ции разностных краевых задач L2m h и Lh натурального числа m. ki , i = 1, 2, . . . , 2. При четном k вклады функций x, y в величины невязок δf1h n - 1, сравнимы между собой; при нечетном k обусловленные функцией x составляющие невязок как имеющие более низкий порядок играют главенствующую роль по сравнению с составляющими, обусловленными функцией y. Практически аналогичной оказалась ситуация с невязками ki с тем лишь отличием, что при нечетном k главенствующую роль δf2h во вклады в величины невязок играет функция y. Декларация о финансовых и других взаимоотношениях. Исследование не имело спонсорской поддержки. Автор несет полную ответственность за предоставление окончательной версии рукописи в печать. Окончательная версия рукописи была одобрена автором. Автор не получал гонорар за статью.

Об авторах

Владимир Николаевич Маклаков

Самарский государственный технический университет

Email: makvo63@yandex.ru
(к.ф.-м.н., доц.; makvo63@yandex.ru), доцент, каф. высшей математики и прикладной информатики Россия, 443100, Самара, ул. Молодогвардейская, 244

Список литературы

Дополнительные файлы