Estimation of the order of the matrix method approximation of numerical integration of boundary-value problems for the second order inhomogeneous linear ordinary differential equations

Full Text

Классический метод численного интегрирования краевых задач для линейных неоднородных обыкновенных дифференциальных уравнений второго порядка (ОДУ2) u + p(x)u + q(x)u = f (x), u0 = u0 , un = un , (1) где u(x) - искомое точное решение; p(x), q(x), f (x) - заданные функции, дифференцируемые нужное число раз; u0 , un - заданные числа, использующий конечные разности, имеет аппроксимацию второго порядка [1-6]. Использующие конечные разности методы численного интегрирования краевых задач для дифференциальных уравнений в частных производных также имеют аппроксимацию второго порядка [4-8]. Последнее обусловлено тем, что при аппроксимации производных конечными разностями в разложении искомой функции в ряд Тейлора удерживалось всего три члена разложения. В работе [9] предложен метод, использующий обратные матрицы и позволяющий увеличить число членов (до произвольного натурального) в разложении искомой функции в ряд Тейлора при численном интегрировании первой краевой задачи для линейных неоднородных ОДУ2 с переменными коэффициентами, при этом аппроксимация производных конечными разностями не использовалась. Обобщение метода для смешанной краевой задачи для линейных неоднородных и для нелинейных ОДУ2 дано в [10, 11], соответственно. Однако оценка порядка аппроксимации метода интегрирования дана не была, исследование которой для первой, второй и третьей краевых задач для линейного неоднородного ОДУ2 поставим целью настоящей работы. Далее будем придерживаться принятых в [4] обозначений: 1) D - область интегрирования, ограниченная отрезком [a, b], Dh - узлы сетки, определяемые значениями xi = x0 + ih, i = 1, 2, . . . , n, x0 = a, xn = b, h = (b - a) /n, n + 1 - число узлов сетки; 2) u(x) - непрерывная функция, являющаяся точным решением краевой задачи (1); 3) [u]h - сеточная функция, совпадающая с точным решением в узлах сетки Dh ; 4) u(h) - искомая сеточная функция. Для краткости примем для любой функции обозначение ϕ(xi ) = ϕi , где xi - узел сетки Dh . В дальнейшем опустим индекс h в наименованиях сеточных функций [u]h , u(h) и будем оговаривать особо случаи, в которых будет использоваться непрерывная функция u(x), являющаяся точным решением, при сохранении обозначений u(xi ) = ui для неё в узлах сетки. 144 Оценка порядка аппроксимации матричного метода численного интегрирования . . . Согласно изложенному в [9] методу численного интегрирования первой краевой задачи (1) для каждого внутреннего узла xi , i = 1, 2, . . . , n - 1 сетки Dh составляется система линейных алгебраический уравнений (СЛАУ), (k) состоящая из k + 1 уравнений относительно неизвестных ui , ui , ui , . . . , ui . k , Pk Первые два уравнения системы есть многочлены Тейлора Pi-1 i+1 степени k, полученные из соответствующих рядов Тейлора: ui-1 = ui - hui + h2 h3 ui - u + ··· = 2! 3! i k (-1)m m=0 ∞ (-1)m + m=k+1 ui+1 = ui +hui + h2 h3 ui + ui +· · · = 2! 3! k m=0 hm (m) u + m! i hm (m) u + m! i hm (m) k k u = Pi-1 + Ri-1 , m! i ∞ m=k+1 hm (m) k k u = Pi+1 +Ri+1 , m! i где k Ri-1 = hk+1 (k+1) u (ξi ) = O(hk+1 ), (k + 1)! ξi ∈ (xi-1 , xi ) - остаточный член. Здесь и далее верхний индекс означает степень многочлена Тейлора, если речь не идёт о показателях алгебраических степеней, степенях производных и символов обратных матриц. Оставшиеся k - 1 уравнений системы есть равенства, полученные дифференцированием по x обеих частей уравнения (1), т.е. q(x)u + p(x)u + u (r) = f (x)(r) , где r = 0, 1, . . . , k-2, которые записаны в узле xi . Преобразования с использованием обратных матриц полученных систем для каждого узла xi позволяют составить систему разностных уравнений для трёхточечного шаблона xi-1 , xi , xi+1 , i = 1, 2, . . . , n - 1, решение которой даёт приближённое искомое решение. В [9] показано, что увеличение степени используемых многочленов Тейлора приводит к уменьшению погрешности между точным и приближённым решениями задачи в узлах сетки. В дальнейшем систему разностных уравнений при фиксированном k буk дем обозначать, по аналогии с [4], в символьной форме как Lk u = fh , причём h k наряду с Lk u = fh , для краткости эту задачу будем обозначать также как h k. Lh 1. Некоторые формальные преобразования и оценки. Для функции u выполним формально следующие преобразования для тройки узлов xi-1 , xi , xi+1 , i = 1, 2, . . . , n - 1, сетки Dh при фиксированном k. Имеем точные равенства [ui ] - h[ui ] + h2 h3 hk (k) k [ui ] - [ui ] + · · · + (-1)k [ui ] = [ui-1 ] - Ri-1 , 2! 3! k! (2) 145 М а к л а к о в В. Н. [ui ] - h[ui ] + h2 hk-1 (k) k-1 [ui ] + · · · + (-1)k-1 [u ] = [ui-1 ] - Ri-1 , 2! (k - 1)! i (3) и h3 hk (k) h2 k [ui ] + [ui ] + · · · + [u ] = [ui+1 ] - Ri+1 , 2! 3! k! i h2 hk-1 (k) k-1 [ui ] + h[ui ] + [ui ] + · · · + [ui ] = [ui+1 ] - Ri+1 . 2! (k - 1)! [ui ] + h[ui ] + (4) (5) Умножая обе части равенств (2), (3) на некоторые числа α0 = 0, β0 = 0 соответственно и складывая, получим: h3 h2 h2 - β0 h [ui ] + -α0 + β0 [ui ]+ 2! 3! 2! hk hk-1 (k) k α0 - β0 [ui ] = α0 [ui-1 ] - Ri-1 + k! (k - 1)! α0 [ui ] + (-α0 h + β0 )[ui ] + α0 + · · · + (-1)k k-1 k-1 k + β0 [ui-1 ] - Ri-1 = [zi-1 ] - α0 Ri-1 - β0 Ri-1 , (6) где введено обозначение [zi-1 ] = α0 [ui-1 ] + β0 [ui-1 ]. (7) Выполняя аналогичные действия с равенствами (4), (5) и числами α1 = 0, β1 = 0, найдём h2 h3 h2 + β1 h [ui ] + α1 + β1 [ui ]+ 2! 3! 2! hk hk-1 (k) k + · · · + α1 + β1 [ui ] = α1 [ui+1 ] - Ri+1 + k! (k - 1)! α1 [ui ] + (α1 h + β1 )[ui ] + α1 k-1 k-1 k + β1 [ui+1 ] - Ri+1 = [zi+1 ] - α1 Ri+1 - β1 Ri+1 , (8) где [zi+1 ] = α1 [ui+1 ] + β1 [ui+1 ]. (9) Для каждого узла xi , i = 1, 2, . . . , n-1, составим СЛАУ, в которую внесём соотношения (6), (8) и равенства q[u] + p[u ] + [u ] 146 (r) = f (r) , r = 0, 1, . . . , k - 2, Оценка порядка аппроксимации матричного метода численного интегрирования . . . записанные в узле xi . В итоге получим   hk-1 hk (k)   α0 [ui ] + (-α0 h + β0 )[ui ] + · · · + (-1)k α0 - β0 [ui ] =   k! (k - 1)!    k-1 k   = [zi-1 ] - α0 Ri-1 - β0 Ri-1 ,     k k-1   (k)  α [u ] + (α h + β ) [u ] + · · · + α h + β h  1 i [u ] =  1 1 1 1 i k! (k - 1)! i k-1 k  = [zi+1 ] - α1 Ri+1 - β1 Ri+1 ,     qi [ui ] + pi [ui ] + [ui ] = fi ,     q [u ] + q + p [u ] + p [u ] + [u ] = f ,  i i  i i i i i i i   . . .     (k-2) (k) (k-2) q [u ] + · · · + [u ] = f . i i i (10) i В матричной форме система (10) имеет вид Ak [Vik ] = [Gk ] i i в обозначениях  -α0 h + β0  α0     α1 α1 h + β1  Ak =  i pi  qi  qi + pi  qi  ...  ... (k-2) qi ... и   [ui ]  [ui ]     [ui ]  k  [Vi ] =  [u ]  ,   i   ...  (k) [ui ] ... ... ... ... ... ... hk-1 hk (-1)k α0 - β0 k! (k - 1)! hk-1 hk α1 + β1 k! (k - 1)! 0 0 ... 1              (11) k-1  k [zi-1 ] - α0 Ri-1 - β0 Ri-1 k  [zi+1 ] - α1 Ri+1 - β1 Rk-1  i+1     fi k . [Gi ] =    fi     ...  (k-2) fi -1 В предположении существования обратной матрицы Ak от матрицы i k найдём Ak -1 [Gk ] = [V k ]. Выпишем первое уравнение последнего матAi i i i ричного равенства k-1 k-1 k k aki [zi-1 ] - α0 Ri-1 - β0 Ri-1 + aki [zi+1 ] - α1 Ri+1 - β1 Ri+1 + 11 12 k+1 (m-3) + aki fi + 13 aki fi 1m = [ui ], m=4 147 М а к л а к о в В. Н. где aki - соответствующие элементы обратной матрицы Ak 1m i или после преобразований [ui ] aki 1 aki - 11 [zi-1 ] + ki - 12 [zi+1 ] = fi + ki ki ki a13 a13 a13 a13 + k -aki (α0 Ri-1 + 11 -1 для узла xi , k+1 (m-3) aki fi 1m m=4 k-1 β0 Ri-1 ) - aki 13 + k-1 k aki (α1 Ri+1 + β1 Ri+1 ) 12 . (12) В дальнейшем матрицы Ak будем называть локальными матрицами, разi мерность которых равна (k + 1). Найдём следующие предварительные оценки: ∞ k Ri+1 - k Ri-1 = m=k+1 hm (m) u - m! i ∞ (-1)m m=k+1 hm (m) u = m! i ∞ = m=k+1 hm (m) u (1 - (-1)m ) (13) m! i и, аналогично, ∞ k k Ri+1 + Ri-1 = m=k+1 hm (m) u + m! i ∞ (-1)m m=k+1 hm (m) u = m! i ∞ = m=k+1 hm (m) u (1 + (-1)m ) . (14) m! i Из равенств (13), (14) для чётного k имеем k k Ri+1 - Ri-1 = k k Ri+1 + Ri-1 = hk+1 (k+1) u (1 + 1) + (k + 1)! i hk+2 (k+2) + u (1 - 1) + · · · = O(hk+1 ), (15) (k + 2)! i hk+1 (k+1) ui (1 - 1) + (k + 1)! hk+2 (k+2) + u (1 + 1) + · · · = O(hk+2 ) (16) (k + 2)! i и для нечётного k - k k Ri+1 - Ri-1 = 148 hk+1 (k+1) ui (1 - 1) + (k + 1)! Оценка порядка аппроксимации матричного метода численного интегрирования . . . + k k Ri+1 + Ri-1 = hk+2 (k+2) u (1 + 1) + · · · = O(hk+2 ), (17) (k + 2)! i hk+1 (k+1) u (1 + 1) + (k + 1)! i hk+2 (k+2) + ui (1 - 1) + · · · = O(hk+1 ). (18) (k + 2)! 2. Порядок аппроксимации первой краевой задачи. Сеточная функция ui , i = 0, 1, . . . , n, являющаяся решением некоторой разностной краевой задачи, при подстановке в уравнения этой разностной краевой задачи обратит их в верные равенства. В [4] показано, что подстановка в уравнения задачи сеточной функции [ui ], отличающейся от ui , приведёт к некоторому отличию k от верных равенств. Эти отличия и характеризует невязка δfh [4]. Иными k словами, подстановка [u] в Lk u = fh приведёт к h k k Lk [u] = fh + δfh . h (19) В соответствии с [4] в качестве оценки величины невязки примем норму k k k k δfh = max |δfh0 |, |δfhn |, max |δfhi | , i = 1, 2, . . . , n - 1, где первые две компоненты характеризуют меру отличий в граничных узлах сетки Dh , оставшиеся - во внутренних узлах. Согласно [4, 5] разностная краевая задача аппроксимирует дифференциk альную краевую задачу на точном решении u, если δfh → 0 при h → 0. Если k , где C > 0, k > 0 - некоторые k Ch при этом имеет место неравенство δfh постоянные, не зависящие от h, то говорят, что имеет место аппроксимация порядка k относительно величины h. k Пусть Lk u = fh аппроксимирует первую дифференциальную краевую заh дачу (1). В этом случае положим в равенствах (7), (9), (12) α0 = 1, β0 = 0 и α1 = 1, β1 = 0 для каждого i = 1, 2, . . . , n - 1. При i = 1 равенство (7) обратится в краевое условие [z0 ] = u0 = u0 , равенство (9) - в [z2 ] = [u2 ], а уравнение (12) примет вид ak1 [u1 ] ak1 1 11 - k1 u0 + k1 - 12 [u2 ] = f1 + k1 k1 a13 a13 a13 a13 k+1 (m-3) ak1 f1 1m m=4 + k k -ak1 R0 - ak1 R2 11 12 . (20) ak1 13 Аналогично при i = n - 1 получено краевое условие [zn ] = un = un , равенство [zn-2 ] = [un-2 ] и уравнение - ak, n-1 11 [un-1 ] ak, n-1 [un-2 ] + k, n-1 - 12n-1 un = ak,n-1 a13 ak, 13 13 149 М а к л а к о в В. Н. = fn-1 + k+1 1 (m-3) ak, n-1 13 ak, n-1 fn-1 1m + k k -ak, n-1 Rn-2 - ak, n-1 Rn 11 12 ak, n-1 13 m=4 . (21) Для каждого i = 2, 3, . . . , n - 2 уравнение (12) запишем так: - [ui ] aki aki 11 [ui-1 ] + ki - 12 [ui+1 ] = aki a13 aki 13 13 1 = fi + ki a13 k+1 (m-3) aki fi 1m + m=4 k k -aki Ri-1 - aki Ri+1 11 12 , aki 13 i = 2, 3, . . . , n - 2. (22) Составим СЛАУ из уравнений (20)-(22). Отбрасывание последних дробей в уравнениях системы, что равносильно переходу от точного решения [ui ] к искомому приближённому ui , приведёт эту систему к системе разностных уравнений, соответствующей первой дифференциальной краевой задаче, по крайней мере при k равном трём и пяти, как и было получено в [9]. Следовательно, в соответствии с (19), последняя дробь в уравнениях (20)-(22), да и в (12), характеризует величину невязки в узлах xi , i = 1, 2, . . . , n - 1, т.е. для рассматриваемой задачи имеем k k -aki Ri-1 - aki Ri+1 11 12 . (23) aki 13 Для первой краевой задачи величина невязки (23) в граничных узлах сетки обращается в нуль в силу того, что уравнения (20), (21) содержат известные значения u0 , un искомой функции в этих узлах. Непосредственными вычислениями убедимся в справедливости оценки k δfhi = k 2 M11 ≈ M11 , (24) k где M11 - алгебраическое дополнение элемента tk транспонированной ло11 кальной матрицы Ak . Действительно, для произвольного натурального числа i k 3 имеем, пренебрегая старшими степенями, k M11 = 150 (-1)k+1 hk k! h h2 2! ... = hk-1 (k - 1)! hk k! pi qi + pi 1 pi ... ... 0 0 pi qi + p i . . . si ui vi wi 1 pi ... ... ... ... ... ... 0 0 ... 1 pi 0 0 ... 0 1 . . . si ui . . . vi wi + ... ... ... . . . 1 pi h h2 2! ... hk-1 (k - 1)! pi = qi + p i . . . si 1 pi ... vi ... ... ... ... = 0 0 ... 1 Оценка порядка аппроксимации матричного метода численного интегрирования . . . k-1 k-1 k-1 cm hm = bk hk + M11 ≈ M11 , = bk hk + m=1 где bk , cm - коэффициенты, не зависящие от h; si , ui , vi , wi - некоторые функции от qi , pi и их производных. Повторное использование несколько раз последней формулы приводит к (24). Аналогично доказывается справедливость оценки k 2 M12 ≈ M12 . (25) k Рассмотрим M13 . Имеем k M13 = + -h h2 2! ... h qi + pi . . . si 2 h pi . . . vi 2! ... ... ... ... hk-1 0 ... 1 (k - 1)! hk 0 ... 0 k! qi + pi . . . si ui (-1)k-1 hk-1 (k - 1)! (-1)k hk k! h h2 -hk 2! = ... k! hk-1 (k - 1)! (-1)k hk k! k-1 = hk m=1 ... ... ... ... + 0 ... 1 pi qi + pi ... si ui pi ... vi wi ... ... ... ... + 0 ... 1 pi h h2 2! ... hk-1 (k - 1)! qi + pi ... si pi ... vi ... ... ... = 0 ... 1 2k-1 k-1 k-1 bm hm + M13 ≈ M13 , cm hm = m=3 1 wi 2k-3 bm hm + = pi vi (-1)k-1 hk-1 (-1)k-1 hk-1 (k - 1)! ... ... (k - 1)! + wi pi -h h2 2! ... -h h2 2! ... ui m=k+1 где bm , cm - коэффициенты, не зависящие от h. После повторного использования несколько раз последней формулы получим k 2 M13 ≈ M13 . (26) 151 М а к л а к о в В. Н. На основании соотношений (24)-(26) для любого k венств Mk ak Mk ak 12 11 - k = - 11 , - k = - 12 k k a13 M13 a13 M13 следуют оценки - M2 ak 2 - hpi 11 ≈ - 11 = , 2 k 2h2 M13 a13 - 3 из очевидных ра- M2 2 + hpi ak 12 ≈ - 12 = . 2 k 2h2 M13 a13 Тогда величину невязки (23) на точном решении [u] во внутренних узлах xi сетки Dh задачи Lk запишем так: h k δfhi = k k k k -ak Ri-1 - ak Ri+1 (2 - hpi ) Ri-1 + (2 + hpi ) Ri+1 11 12 = = 2h2 ak 13 = k k k k Ri+1 + Ri-1 hpi Ri+1 - Ri-1 + h2 2h2 (27) для произвольного k. Для чётного k с учётом (15), (16) из (27) имеем k δfhi = k k k k Ri+1 + Ri-1 hpi Ri+1 - Ri-1 + = h2 2h2 O(hk+2 ) hpi O(hk+1 ) = + = O(hk ) + O(hk ) = O(hk ), (28) h2 2h2 и для нечётного k с учётом (17), (18) - k δfhi = O(hk+1 ) hpi O(hk+2 ) + = O(hk-1 ) + O(hk+1 ) = O(hk-1 ). h2 2h2 (29) k Из равенства (28) следует оценка δfh Chk для чётного k, а из (29) - k k-1 для нечётного k, откуда имеем, что задачи L2m и L2m+1 оценка δfh Ch h h для любого натурального m 1 имеют одинаковый порядок аппроксимации. Следовательно, на практике для уменьшения объёма вычислений следует использовать задачи Lk с чётными значениями k. Действительно, число арифh метических операций только для нахождения обратной матрицы от локальной Ak методом Гаусса вычисляется по формуле [12] i (k + 1)2 k + 1 8 (k + 1)3 - - - 1. 3 2 6 Результат о выборе чётного значения k ниже будет подтверждён численным экспериментом. 3. Порядок аппроксимации второй и третьей краевых задач. При исследовании второй и третьей краевых задач в [4-6] производные в граничных узлах сетки были заменены конечными разностями первого порядка аппроксимации относительно h, в силу чего СЛАУ для вычисления ui стала содержать 152 Оценка порядка аппроксимации матричного метода численного интегрирования . . . n + 1 уравнений с n + 1 неизвестными (по сравнению с первой краевой задачей в систему были добавлены два уравнения, содержащие неизвестные u0 , un ). В силу того, что СЛАУ для вычисления ui , составленная из уравнений (20)-(22), записанных без остаточных членов, содержит n-1 уравнений с n- k - 1 неизвестными, то не представляется возможным вычислить невязки δfh0 k с использованием последней дроби в (12). Поэтому далее в качестве и δfhn k k нормы примем δfh = max |δfhi | , i = 1, 2, . . . , n - 1. В узлах сетки Dh с номерами i = 2, 3, . . . , n - 2 второй и третьей краевых задач, граничные условия которых имеют соответственно вид u0 = u0 , un = un и α0 u0 +β0 u0 = z0 , α1 un +β1 un = zn , где u0 , un , z0 , zn , α0 , α1 , β0 , β1 - заданные числа, ситуация полностью совпадает с изложенной выше первой краевой задачей, для которой α0 = 1, β0 = 0 и α1 = 1, β1 = 0. Оценим порядок аппроксимации третьей краевой задачи в узле x1 , в котором α0 = 0, β0 = 0 и α1 = 1, β1 = 0. Из последней дроби равенства (12) следует оценка невязки k δfh1 = k-1 k k -ak1 (α0 R0 + β0 R0 ) - ak1 R2 12 11 . ak1 13 (30) Локальная матрица Ak рассматриваемой задачи отличается только вто1 рой строкой, составленной из коэффициентов разложения (4), от локальной k матрицы (11) при i = 1, поэтому оценка (24) для M11 сохранит форму, а для k k M12 , M13 найдём: k-1 k-1 2 k M12 = bk hk + ck hk-1 + M12 ≈ M12 ≈ M12 , 2k-1 k M13 2k-2 k-1 k-1 2 cm hm + M13 ≈ M13 ≈ M13 , m = bm h + m=k+1 m=k+1 что совпадает с (25) и (26) соответственно. На основании соотношений (24)-(26) для любого k венств ak1 Mk ak1 Mk 11 - k1 = - 11 , - 12 = - 12 k k a13 M13 ak1 M13 13 следуют оценки: - - 3 из очевидных раak1 M2 2 - hp1 2 11 ≈ - 11 = ≈ , 2 2 k1 -3β0 h + 2α0 h -3β0 h M13 a13 ak1 M2 -2β0 - 2h (β0 p1 - α0 ) + α0 h2 p1 2 12 ≈ - 12 = ≈ 2. 2 2 + 2α h3 k1 -3β0 h 3h M13 a13 0 Тогда из соотношения (30) имеем k δfh1 = 2 α0 O(hk+1 ) + β0 O hk -3β0 h + 2 O(hk+1 ) = 3h2 = O(hk-1 ) + O(hk-1 ) = O(hk-1 ), 153 М а к л а к о в В. Н. т.е. порядок аппроксимации в узле x1 , следовательно, и всей разностной краевой задачи Lk оказался на единицу меньше степени многочлена Тейлора k. h Точно такой же вывод о степени аппроксимации второй разностной краевой задачи следует из равенства (12), записанного в узле x1 , в котором α0 = 0, β0 = 1 и α1 = 1, β1 = 0: k δfh1 = k-1 k -ak1 R0 - ak1 R2 11 12 . ak1 13 (31) Ситуация в узле xn-1 третьей (α0 = 1, β0 = 0, и α1 = 0, β1 = 0) и второй (α0 = 1, β0 = 0, и α1 = 0, β1 = 1) краевых задач полностью аналогична изложенной. Действительно, в этом случае локальная матрица Ak задачи n-1 будет отличаться только первой строкой, составленной из коэффициентов разложения (2), от локальной матрицы (11) при i = n - 1. Невязка (31) показывает, что порядок аппроксимации второй разностной краевой задачи L2 равен единице в узлах x1 , xn-1 и, как следует из (28) для h чётного k = 2, равен двум в узлах i = 2, 3, . . . , n-2. Повысим порядок аппроксимации на единицу в узлах x1 , xn-1 , оставаясь в рамке второй разностной краевой задачи L2 . h Запишем систему (10) в узле x1 при α0 = 0, β0 = 1 и α1 = 1, β1 = 0, ограничиваясь степенями производных не старше третьей:  2  [u ] - h[u ] + h [u ] = u - R2 ,  1  0 0 1  2! 1  h2 h3 3 [u1 ] + [u ] = [u2 ] - R2 ,  [u1 ] + h[u1 ] +  2! 3! 1    q1 [u1 ] + p1 [u1 ] + [u1 ] = f1 , в которую подставим значение [u1 ], найденное дифференцированием обеих частей равенства: [u ] = f - p[u ] - q[u]. (32) Получим  2 2  h  - q1 [u1 ] + 1 - h q1 + p1 [u1 ] - h +   2! 2!    h3 h3  1- q1 [u1 ] + h - q1 + p1 [u1 ] + 3! 3!         q1 [u1 ] + p1 [u1 ] + [u1 ] = f1 . h2 h2 2 p1 [u1 ] = u0 - f - R0 , 2! 2 2! 1 h h3 - p1 [u1 ] = 2! 3! h3 3 = [u2 ] - f - R2 , 3! 1 В матричной форме последняя система уравнений имеет вид A2 [V12 ] = [G2 ] 1 1 154 Оценка порядка аппроксимации матричного метода численного интегрирования . . . в обозначениях h2 - q1  2!  2 3 A1 =  1 - h q  3! 1 q1  h2 1- q1 + p 1 2! h3 h- q1 + p1 3! p1 h2 -h - p1 2! h2 h3 - p1 2! 3! 1    ,    h2 2  u0 - 2! f1 - R0    2 h3 [G1 ] =  3 . f1 - R2   [u2 ] - 3! f1    [u1 ] [V12 ] =  [u1 ] , [u1 ] -1 В предположении существования обратной матрицы A2 от локаль1 2 найдём A2 -1 [G2 ] = [V 2 ]. Выпишем первое уравнение последнего ной A1 1 1 1 матричного равенства a2 u0 - 11 h3 h2 2 3 f1 - R0 + a2 [u2 ] - f - R2 + a3 f1 = [u1 ], 12 13 2! 3! 1 или - 2 3 a 2 h3 -a2 R0 - a2 R2 a2 [u1 ] a2 a2 h2 12 11 12 11 11 f1 - 2 f1 + u0 + 2 - 11 [u2 ] = f1 - 2 , a2 a13 a2 a13 2! a13 3! a2 13 13 13 где последняя дробь характеризует величину невязки. Для рассматриваемого случая имеем следующие оценки коэффициентов: - -12 + 6hp1 - 2h2 p2 - q1 - p1 a2 -2 1 11 = ≈ , 2 2 p - 5h3 (q + p ) 18h + 4h 1 3h a13 1 1 - 12 + 12hp1 + 6h2 p2 - q1 - p1 a2 2 1 12 = ≈ 2, 2 2 + 4h3 p - 5h4 (q + p ) 18h 3h a13 1 1 1 тогда -2 O h3 2 O(h4 ) + = O(h2 ) + O(h2 ) = O(h2 ), 3h 3h2 т.е. порядок аппроксимации в узле x1 и, аналогично, в xn-1 , следовательно, и всей разностной краевой задачи L2 теперь оказался равным двум. h Отметим, что такой же (второй) порядок аппроксимации имеет вторая разностная краевая задача L3 . Действительно, из (29) для нечётного k = 3 h следует второй порядок аппроксимации в узлах xi , i = 2, 3, . . . , n - 2. Точно такой же (второй) порядок даёт невязка (31) при k = 3 в узлах x1 , xn-1 . Таким образом, из вышеизложенного следует, что увеличение степени многочлена Тейлора k и использование операции дифференцирования обеих частей равенства (32) k - 1 раз позволяют аналогичным образом увеличить порядок аппроксимации в узлах x1 , xn-1 , следовательно, и всей задачи до произвольного натурального числа. 2 δfh1 = 155 М а к л а к о в В. Н. Использование в системе (10) различных комбинаций значений α0 , β0 и α1 , β1 для узлов x1 и xn-1 соответственно даёт возможность рассмотреть краевые задачи со смешанными краевыми условиями [10]. 4. Оценка погрешности. Для задачи Lk были приняты следующие две h нормы для погрешности в узлах xi , i = 1, 2, . . . , n - 1 сетки Dh : n-1 2 n-1 1) Dk = 100% - суммарная оценка отноi=1 ([ui ] - ui ) i=1 [ui ] сительной погрешности, которую можно трактовать как некий аналог коэффициента вариации в статистике, который характеризует меру разброса в процентах [13]; величина Dk отличается от коэффициента вариации тем, что стандартное отклонение заменено корнем квадратным из остаточной дисперсии [13]. 2) E k = max [ui ]-ui , i = 1, 2, . . . , n-1 - максимальная оценка абсолютной погрешности. В качестве примера использовано ОДУ2 u - 2 2 u + 2 = x cos x, x x (33) имеющее аналитическое решение u = C1 x + C2 x2 - x cos x, и исследованы первая и третья краевые задачи с краевыми условиями u(5) = 22.32, u(13) = 23.95 (34) и u(5) + 3u (5) = 17.60, 2u(13) + 2u (13) = 56.02 (35) соответственно. Было принято n = 20, h = 0.4. Результаты численного эксперимента для первой краевой задачи (33), (34) приведены в табл. 1, для третьей - в табл. 2. Расчёты для третьей краевой задачи (33), (35) выполнены без использования метода повышения порядка аппроксимации. Анализ данных табл. 1 свидетельствует, что суммарная относительная Dk 2m+1 и максимальная абсолютная E k погрешности задач L2m и Lh , имеющих h одинаковый порядок аппроксимации, различаются незначительно для любого натурального m ∈ [1, 4], тогда как указанная особенность в данных табл. 2 отсутствует, чего и следовало ожидать. Отметим, что погрешности для первой краевой задачи L2m и третьей краевой задачи L2m+1 сравнимы между h h собой, что подтверждает вывод о том, что порядок аппроксимации второй и Таблица 1 Значения погрешностей для первой краевой задачи (33), (34) [The accuracies for the first boundary value problem (33), (34)] k 2 3 7.29 · 10 1.94 · 10-1 4.84 · 10 1.30 · 10-3 4.18 · 10-4 1.04 · 10-3 k 6 7 8 9 D ,% Ek 156 -6 1.96 · 10 5.37 · 10-6 -2 5 8.11 · 10 2.23 · 10-1 k -2 4 D ,% Ek k -6 1.51 · 10 4.10 · 10-6 -4 -7 5.88 · 10 1.21 · 10-6 5.90 · 10-7 1.21 · 10-6 Оценка порядка аппроксимации матричного метода численного интегрирования . . . Таблица 2 Значения погрешностей для третьей краевой задачи (33), (35) [The accuracies for the third boundary value problem (33), (35)] k k D ,% Ek k k D ,% Ek 2 3 -1 4 5 6 -5 6.88 · 10 1.15 · 10-4 1.08 · 10 2.79 · 10-1 1.06 · 10 1.79 · 10-2 4.22 · 10-4 1.18 · 10-3 7 6.42 · 10 1.16 -1 8 9 -6 1.24 · 10 2.80 · 10-6 -2 -7 2.98 · 10 7.10 · 10-7 3.04 · 10-7 6.80 · 10-7 третьей разностных краевых задач на единицу меньше степени используемого многочлена Тейлора. Выводы. 1. Теоретически выявлены закономерности между порядком аппроксимации метода и степенью используемого многочлена Тейлора для краевых задач с граничными условиями различных родов. Установлено, что порядок 2m+1 совпадааппроксимации для первых разностных краевых задач L2m и Lh h ет и равен 2m для любого натурального m, тогда как для второй и третьей краевых задач указанная особенность отсутствует, причем порядок аппроксимации задач прямо пропорционален степени используемого многочлена Тейлора и меньше него на единицу. 2. Предложен способ повышения порядка аппроксимации метода до произвольного натурального числа для разностных краевых задач с граничными условиями второго и третьего рода.

About the authors

Vladimir N Maklakov

Samara State Technical University

Email: makvo63@yandex.ru
(Cand. Phys. & Math. Sci.; makvo63@yandex.ru), Associate Professor, Dept. of Higher Mathematics and Applied Informatics 244, Molodogvardeyskaya st., Samara, 443100, Russian Federation

References

Supplementary files