Вестник КРАУНЦ. Физ.-мат. науки. 2021. Т. 34. №1. C. 88-104. ISSN 2079-6641
УДК 519.644 Научная статья
О приближении определенных интегралов составными квадратурными формулами с использованием производных
В. В. Шустов
Государственный научно-исследовательский институт авиационных систем, 125319, г. Москва, ул. Викторенко, дом 7 E-mail: [email protected]
Рассмотрена задача вычисления определенного интеграла функции, для которой известны значения ее самой и набора производных до заданного порядка в точках отрезка интегрирования. Построены составные квадратурные формулы, которые используют значения функции и ее производных до m-го порядка включительно. Получено представление остаточного члена, выраженное через производную соответствующего порядка и число узловых точек. Приведены примеры интегрирования заданных функций с исследованием погрешности и ее оценки. Дано сравнение с известными численными методами и формулой Эйлера-Маклорена, которое показало повышенную точность и лучшую сходимость метода двухточечного интегрирования.
Ключевые слова: квадратурные формулы с использованием производных, составная квадратурная формула, остаточный член интегрирования, оценка погрешности приближения, формула Эйлера-Маклорена
DOI: 10.26117/2079-6641-2021-34-1-88-104
Поступила в редакцию: 29.12.2020 В окончательном варианте: 21.03.2021
Для цитирования. Шустов В. В. О приближении определенных интегралов составными квадратурными формулами с использованием производных // Вестник КРАУНЦ. Физ.-мат. науки. 2021. Т. 34. № 1. C. 88-104. DOI: 10.26117/2079-6641-2021-34-1-88-104
Контент публикуется на условиях лицензии Creative Commons Attribution 4.0 International (https://creativecommons.org/licenses/by/4.0/deed.ru)
© Шустов В. В., 2021
Введение
Известными методами вычисления определенных интегралов являются методы трапеций, Симпсона, Гаусса, Чебышева и другие, изложенные в учебниках общего плана [1]-[7], более специализированных изданиях [8]-[9] и других работах. В этих методах для приближения определенных интегралов используются только значения функции, заданных в узловых точках отрезка интегрирования. В тоже время для функций, обладающих определенной степенью гладкости, представляет интерес и является важным в теоретическом и практическом отношении использование
Финансирование. Исследование выполнялось без финансовой поддержки фондов.
также значений ее производных в этих точках. Известным методом интегрирования, в котором используется значения производных, является метод, основанный на применении формулы Эйлера-Маклорена. Однако в этом методе, как отмечено в [10, с. 544] получается ряд, который, вообще говоря, расходится.
В работах [11] и [16] предложен метод интегрирования, основанный на использовании двухточечных многочленов Эрмита в варианте, когда производные заданы только на концах отрезка интегрирования.
В рамках развития предложенного метода в настоящей работе рассматривается случай, когда производные заданы не только на концах отрезка интегрирования, но и в узловых точках, расположенных внутри этого отрезка. В этом случае используется подход, согласно которому весь отрезок интегрирования разбивается узловыми точками на элементарные отрезки, на границах которых заданы значения функции и ее производных до некоторого порядка, а также используется аддитивность операции интегрирования относительно промежутка интегрирования.
Необходимость получения и использования составных квадратурных формул обусловлена тем, что в ряде случаев трудно или невозможно получить данные о производных высокого порядка, или сходимость процесса интегрирования не является удовлетворительной.
Задача интегрирования функции со многими узловыми точками сводится к аналогичной задаче только с двумя точками, расположенными на концах отрезка интегрирования. Эта задача решена и представлена в работах [16] и [17], в которых для построения квадратурных формул используется двухточечные интерполяционные многочлены Эрмита [12]. При этом рассматривается общий случай, когда узловые точки могут быть расположены на разных расстояниях друг от друга.
1. Постановка и решение задачи
Пусть функция f (х) определена на отрезке [х0,хп] и является дифференцируемой на этом отрезке достаточное число раз. Пусть также в п + 1 узловых точках хг-, I = 0,1,..., п, таких что Хо < Х1 < Х2 ■ ■ ■ < Хп
этого отрезка заданы значения функции f(x) и ее производных до порядка т включительно:
/Л(Х1)= ] = 0,1,...,т, I = 0,1,...,п. (1)
Из условия существования производных следует, что для функции f(x) существует интеграл I, определенный формулой
Г Хп
I = f (Х)йХ. (2)
Jxо
Необходимо построить квадратурную формулу для интеграла 1т, который приближает интеграл I в соответствии с формулой
ГХп
/ f (x)dx = 1т + Гт, (3)
■'Х0
и использует условия (1), а также получить и оценить ее остаточный член приближения гт.
В силу аддитивности операции интегрирования, в соответствии которой интеграл по отрезку [хо,хп] равен сумме интегралов по составляющим его локальным отрезкам [х,,хг+1],г = 0,1,...,п — 1, имеет место формула
гх* П-11 Г Х!+1
/ / (х)^х = £ / / (х)^х. (4)
./Хо г=0
Согласно результатам, представленным в [11, с. 85] и [16], интеграл по локальному отрезку [х,,х,+1 ] при условиях (1) в наших обозначениях может быть представлен в виде
№+1
ш
/(x)dX = £ Dmhj+1[/i(j) + (-1) j/i+1] + ri, (5)
j=0
с+1
п; = ст+1
т (7 + 1)!с2++2 ,
= (—1)Щ+1 /(2т+2)(п.)^2т+3 (т + 1)!(т + 1)! (6)
Г = (2т + 2)! (2т + 3)! , (6)
й, = хг-+1 — х, и точки Пг ^ (х,,хг-+1). Здесь ст — биноминальный коэффициент. Формулу (4) для интеграла I с учетом (5) можно записать в виде
Л*П п—1 т
/ / (х)^х = £ (£ П,^1/^—1) />] + гг), (7)
■Ух0 г=0 7=0
а также, раскрывая скобки в (7), в виде
гхп п—1 т , л п—1
/ /(х)^х = £ £ птл/+1 + (— 1)/+1] + £ гг. (8)
■Ух0 г=0 7=0 г=0
Из формул (3) и (7) следует, что интеграл 1т может быть записан в виде
п— 1 т
1т = £ £п>/+/ + (— 1)/Й], г=0 7=0
а остаточный член гт как
п—1
Гт = £ г,-. (9)
=0
Подставляя соотношение (6) в (8), вынося коэффициенты, не зависящие от переменной г перед знаком суммирования по этой переменной и используя (9), получим, что для остаточного члена гт имеет место соотношение
(_т+1 ь п—1
Гт = С(2т)+2)Г £ ^2т+3/(2т+2)(П ,). (Ю)
В соответствии с полученными результатами имеет место следующая теорема.
Теорема 1. Пусть функция f(x) удовлетворяет условиям (1). Тогда для определенного интеграла этой функции имеет место
f (x)dx = Im + Гт, (11)
'Х0
где составная квадратурная формула представления интеграла 1т имеет вид
п-1 т
1т = £ £ От^^^—У]], (12)
1=0 ]=0
С]+1
о] = Ст+1
т (] + 1)С++2,
остаточный член гт есть
У^Е А?"+3 f (Jmt2)(ч-),
Гт = + 2)Г £ НГ+'^2т+2,(П). (13)
коэффициенты двухточечного интегрирования Ьт представляются как
Ь = (т + 1)!(т + 1)! (14)
Ьт = (2т + 3)! , (14)
к = Хг-+1 — Х1 и точки п е (хг-,Хг-+1), I = 0,1,...,п — 1. Здесь С]т — биноминальный коэффициент.
Формулу (12) для интеграла 1т, изменяя порядки суммирования по индексам I и ] и вынося коэффициент, не зависящий от индекса I из-под знака его суммирования, можно представить также в виде
т п— 1
1т = I от £ ^[¿^(—^^]. (15)
]=0 1=0
Частный случай. При равномерном распределении узловых точек, при котором шаг между узловыми точками один и тот же:
к = к, I = 0,1,..., п — 1, формулы интегрирования можно упростить.
Формула (15) для интеграла 1т при этом условии принимает вид:
т п—1
1т = I Оттк]+1 I ^]) + (—1)(16) ]=0 1=0
которую также можно записать как
Im = Е LfW + (-l)Jf+], (17)
j=0 n i=0
где длина Ь = Хп — Х0 отрезка [х0,Хп] связана с шагом к и числом отрезков разбиения п соотношением
к = -. (18)
п
x
Формула (10) для остаточного члена гш запишется как:
Г" = + % /<2"+2,(П'). (19)
Для упрощения суммы, стоящей в правой части (19), применим теорему о среднем для определенного интеграла (например, [13, с. 400]), согласно которой, если функция ¥(х) является непрерывной на отрезке [а,Ь], то существует точка П е (а,Ь), такая, что
f F(x)dx = F(n)(b-a). Ja
Применяя теорему о среднем для интеграла по всему отрезку [хо,хп], можно записать:
/ ¥ (х)Жс = ¥ (п )(хи - хо)= ¥ (п )иН.
Jxо
С другой стороны, применяя теорему о среднем для совокупности интегралов по всем отрезкам [хг-,хг+|] и используя аддитивность операции интегрирования, получим, что
/ ¥(х)^х = % / ¥(х)^х = % ¥(п¿)(хг+' -хг') = Н % ¥(п,■).
■^о г'=о х г=0 г=0
Из сравнения обеих цепочек равенств следует, что имеет место равенство
И—'
% ¥(п■)= «¥(п), '=о
в соответствии с которым выполняется соотношение
И'
£ f (2m+2)(ni)= nf (2m+2)(n).
'=о
С использованием этой формулы выражение (19) для остаточного члена гт запишется в виде:
Г" = (—'СГ3 / ^), (20)
а также с учетом (18) в виде
_ = (—')т+'Ь"¥2"+3 (2"+2)(п) (21) " = (2" + 2)!и2"+2 * (п). (21)
Следствие 1. Пусть производная функции порядка 2" + 2 включительно на отрезке [хо,хИ] ограничена некоторой константой М2Ш+2 > 0, т.е. считаем, что
|/^(х)^ М2"+2, х е (хо,хи).
Тогда для погрешности приближения интеграла 8ш= |г" имеет место
где Ат обозначена оценка погрешности приближения
Am —
bmnh2m+3
7M2m+2,
(2т + 2) Г
которую с использованием (18) также можно записать в виде:
Ьт.Ь2т+3
Am —
(2m + 2) !n2m+2
M2m+2.
(22) Таблица 1
Формулы численного интегрирования
m Формулы для интеграла Im
0 1о — 2 (fo + 2fi + ... + 2 fn-i + fn)
1 Ii — 2(fo + 2fi + ... + 2fn-i + fn) + h2(fo - fn)
2 I2 — h (fo + 2fi + ... + 2 fn-i + fn) + ho (fo - fn) + ih0 (fo + 2fl + ... + 2fn-i + fü)
3 I3 — h (fo + 2 fi +... + 2 fn-i + fn) + % (fo- fn) + f + 2fi +... + 2 fn-i + 0+ 1 h4 ( j/// j/// ) 1 168o (Jo Jn )
4 I4 — h (fo + 2fi... + 2fn-i + fn) + % (fo - fn) + t2 (fo + 2 fi + ... + 2fn- i + 1 i ho8 (fo fn ) I 3o24o (fo + 2fi + ... + 2fn-i + fn )
5 I5 — h (fo + 2 fi... + 2fn-i + fn) + f (f/-fn) + 66 (fo + 2fi + ... + 2f- + f0 + + 792 (fo fn )+ 1584o (fo 1 2fi 1 ". 1 2fn-i 1 fn) 1 66h28o (fo fn)
В табл. 1 приведены формулы интегрирования, которые являются частными случаями общей квадратурной формулы (15) для начальных значений т.
Из формулы (12) и из формул, представленных в табл. 1, видно, что каждая формула является суммой компонент, каждая из которых отвечает определенному значению порядка т используемых производных. При этом для четных т в каждой компоненте суммируются значения производных во всех узловых точках, причем для внутренних точек эти значения удваиваются. Для нечетных же значений т соответствующий компонент использует разность значений производных только в концевых точках, т.к. в силу структуры формулы (12) сумма значений производных нечетного порядка во внутренних точках обнуляется.
Можно отметить, что формула, размещенная в строке со значением т=0, соответствует составной квадратурной формуле метода трапеций. Таблица 2
Выражения для оценки погрешности интегрирования А1т
m i (2m+2)! bm AIm
0 i 2 i 6 M2nh3 i2
1 i 24 i 3o M4nh5 72 o
2 i 72 o i 14o M6nh' ioo8oo
3 i 4o32o i 63o M8nh9 254oi6oo
4 i 36288oo i 2772 Mio nh11 ioo59o336oo
5 i i Mi2nhi3
479ooi6oo 12o12 57537672192oo
6 i i Mi4nhi5
87i7829i2oo 5i48o 448793843o976ooo
В табл. 2 приведены выражения для оценки погрешности интегрирования А/ш и компоненты ее представления для начальных значений параметра ".
Из табл. 2 видно, что оценка погрешности интегрирования при "=0, представленная формулой А/о = ^¡^ , соответствует оценке погрешности метода трапеций (см., например, [4, с. 107]).
2. Результаты численного эксперимента
В качестве примера применения полученных составных квадратурных формул и оценки их точности вычислим значение интеграла от функции f (x) = sinx на отрезке [0, п], т.е.
Í п
I = sin xdx.
Jo
Точное значение этого интеграла легко определяется аналитически и есть I = 2. Используя квадратурные формулы, представленные в табл. 1, и значения производных подынтегральной функции, определяемые соотношением ([13, с.149]):
(sin x)(j) = sin (x + 2 j), (23)
вычислены значения интеграла Im и численной погрешности 8Im для различных значений m и числа отрезков разбиения n.
Численные значения Im и погрешности приближения 8Im для значений m = 0 — 7 и числа отрезков разбиения n представлены в табл. 3.
Таблица 3
Значения интеграла Im и его погрешности 8Im для различных m и n
n m 1 2 4 8 16
0 0.0000000 2.0000000 1.5707963 0.42920367 1.8961189 0.10388110 1.9742316 0.025768398 1.9935703 0.0064296562
1 1.6449341 0.35506593 1.9820298 0.017970156 1.9989273 0.0010727229 1.9999337 0.000066303260 1.9999959 4.1325290x10-6
2 1.9739209 0.026079120 1.9996801 0.00031986290 1.9999953 4.7381119x10-6 1.9999999 7.3078996x10-8 2.0000000 1.1381883x10-9
3 1.9989520 0.0010479748 1.9999968 3.1515877x10-6 2.0000000 1.1616152x10-8 2.0000000 4.4738457x10-11 2.0000000 1.7414686x10-13
4 1.9999734 0.000026583556 2.0000000 1.9722292x10-8 2.0000000 1.8114062x10-11 2.0000000 1.7427003x10-14 2.0000000 1.6955457x10-17
5 1.9999995 4.6462431x10-7 2.0000000 8.5345467x10-11 2.0000000 1.9549848x10-14 2.0000000 4.6992911x10-18 2.0000000 1.1428645x10-21
6 2.0000000 5.9369402x10-9 2.0000000 2.7063220x10-13 2.0000000 1.5470038x10-17 2.0000000 9.2922992x10-22 2.0000000 5.6490617x10-26
7 2.0000000 5.7891324x10-11 2.0000000 6.5591947x10-16 2.0000000 9.3600562x10-21 2.0000000 1.4050592x10-25 2.0000000 2.1352508x10-30
В каждой ячейке таблицы сверху представлены значения интеграла /ш, снизу приведены значения погрешности приближения 8/ш, как модуля разности между
точным и приближенным значениями интегралов, определенной по формуле
S Im =
rx 1 rx 1
/ f (x)dx — Hm(x)dx
Jx0 Jx0
(24)
Из данных, представленных в табл. 3, видно, что как с увеличением числа узловых точек п, так и с увеличением порядка производных т последовательные приближения интеграла 1т сходятся к точному его значению I, однако с разной скоростью. Если для метода трапеций, отвечающего значению т = 0, скорость сходимости является самой малой, то с увеличением порядка производных т скорость сходимости приближений существенно возрастает и обеспечивает необходимую точность уже при малом числе точек разбиения п.
Влияние параметров п и т на скорость сходимости хорошо видно при представлении результатов расчета в графической форме, как это показано на рис. 1 и рис. 2.
Рис. 1. Зависимость оценки погрешности Л(п) для значений т = 0 — 7
Рис. 2. Зависимость погрешности Л(т) для значений п = 1 — 8
На этих рисунках представлена оценка погрешности приближения, определенная формулой (22), с константой М2т+2 = 1, т.е.
д = ЬтЬ2т+3
(2т + 2)!п2т+2'
представленная в параметрической форме Д(п) и Д(т) в виде зависимости ее от п и т, соответственно. Отметим, что из-за большого диапазона изменения для осей ординат используется логарифмическая шкала. Зависимости Д(п) и Д(т) представлены в виде точек, соединенных отрезками прямых.
Из представленных графиков видно, что оценка погрешности монотонно убывает как с возрастанием порядка производных т, так и с увеличением числа узловых точек п, причем скорость убывания погрешности увеличивается с возрастанием т (см. рис. 1).
3. Сравнение результатов, полученных по двухточечной формуле интегрирования и по формуле Эйлера-Маклорена
Интересно сравнить предложенные квадратурные формулы с другими подходами, которые также используют значения функции и ее производных на отрезке интегрирования, в частности, основанными на применении формулы Эйлера-Маклорена.
Представим в этом случае интеграл
Г Хп
I = Л (х)йх Ло
аналогично в виде
1 = Ет + гЕт>
где приближенное значение интеграла Ет выражается формулой Эйлера-Маклорена (например, [14, с. 136]), записанной в наших обозначениях как
Ет = 2Л + 2Л + ... + 2^-1 + Л] + ¿Л-(25) гЕт - остаточный член, который согласно [1, с. 292] представляется как
гЕт = -^Л(2т+2)(Я)'Я 6 (Х0'хп), (26)
где В2] - числа Бернулли, рассмотренные, например, в [1, с. 292]. Числа Бернулли представляется формулой
В'2<» = <-1)"+1 ÜpÜ «<2'")' (27)
где £(•?) — дзета-функция Римана (например, [10, с. 263]), определенная как
^ 1
с (») = £ р ■
к=1 к
и обладающая свойством [1, с. 293], что £ ^) ^ 1 при 5 ^
Формулу Эйлера-Маклорена также можно записать в развернутом виде
Ет = I[/о + 2/1 + ... + 2/я_1 + /п] + ^ /0 - /П]■+ ^ /0" - /п ] + ....+ (28)
, В2тЬ2т й(2т-1) /(2т_1)1 [/о _/п ].
Сравнение формул (16) и (28) показывает, что в отличие от формулы двухточечного интегрирования, в которую входят производные как четного, так и нечетного порядков, формула Эйлера-Маклорена содержит производные только нечетного порядка.
Для погрешности приближения 8Ет = \гЕт\ формулы Эйлера-Маклорена при условии, что производная функции порядка 2т + 2 включительно на отрезке [хо,Х1] ограничена некоторой константой М2т+2 > 0, т.е. считая, что
\/(2т+2)(х)\< М2т+2, х Е (Х0,Х1),
имеет место
8Ет < ДЕт,1
где ДЕт обозначена оценка погрешности приближения интегрирования
ДЕ = |Вт+2|^2т+3 М
ДЕт = (2т + 2)! М2т+2. (29)
С целью сопоставления обоих методов интегрирования проведены расчеты с использованием квадратурной формулы Эйлера-Маклорена для двух тестовых функций.
Численный пример 1
В первом примере проведено вычисление ранее рассмотренного в разделе 2 интеграла
Í п
I = sin xdx, Jo
с использованием формулы (25) и его погрешности, определенной формулой (24).
Результаты расчетов в виде приближенных значений интегралов Em и их численных погрешностей 8Em представлены в в соответствующих столбцах табл. 4. В каждой ячейке этой таблицы вверху приводятся значение интеграла Em , внизу приводится значение его погрешности 8Em. Для сопоставления в этой же таблице приведены значения интеграла Im, полученного по формуле двухточечного интегрирования и его погрешности 8Im, также в соответствуюших столбцах.
Анализ данных, представленных в таблице, показал, что приближения интеграла, полученные по двухточечной формуле и по формуле Эйлера-Маклорена, сходятся к точному значению интеграла как с увеличением m, так и с увеличением п. При этом значения погрешности двухточечной формулы при m > 2 меньше, чем значения погрешности формулы Эйлера-Маклорена.
Таблица 4
Значения интеграла 1т и его погрешности 51т для п = 1,2,4, полученные по двухточечной формуле и по формуле Эйлера-Маклорена
n m
\
новый метод
1
метод
Эйлера-
Маклорена
2
новый метод
2
метод
Эйлера-
Маклорена
4
новый метод
4
метод
Эйлера-
Маклорена
0
0.0000000 2.0000000
0.0000000 2.0000000
1.5707963 0.4292037
1.5707963 0.4292037
1.8961189 0.1038811
1.8961189 0.1038811
1
2
3
4
5
6 7
1.6449341 0.3550659
1.6449341 0.3550659
1.9820298 0.0179702
1.9820298 0.0179702
1.9989273 0.0010727
1.9989273 0.0010727
1.9739209 0.0260791
1.9155149 0.0844851
1.9996801 0.0003199
1.9989411 0.0010589
1.9999953 4.738x10-6
1.9999842 0.0000158
1.9989520 0.0010480
1.9790988 0.0209011
1.9999968 3.152x10-6
1.9999346 6.536х10
-5
2.0000000 1.162x10-8
1.9999998 2.432х10
-7
1.9999734 0.0000266
1.9947875 0.0052125
2.0000000 1.972x10-8
1.9999959 4.073х10
6
2.0000000 1.811x10
11
2.0000000
3.788х10-9
1.9999995 4.646x10-7
1.9986977 0.0013023
2.0000000 8.535x10-1
1.9999997 2.544х10
7
2.0000000 1.955x10-14
2.0000000 5.915х10-11
2.0000000 5.937x10-9
1.9996745 3.2610-4
2.0000000 2.706x10-13
2.0000000 1.590х10-8
2.0000000 1.547x10-1
2.0000000 9.240х10-13
2.0000000 5.79x10
11
1.9999186 8.14 x10
5
2.0000000 6.559x10-16
2.0000000 9.934х10
10
2.0000000 9.360x10
-21
2.0000000 1.444х10
14
í(m)
0 1 2 3 4 i t 7
Рис. 3. Сравнение погрешности по формуле двухточечного интегрирования и Эйлера-Маклорена для функции f (x) = sinx
На рис. 3 данные расчетов представлены в графической форме, которая наглядно показывает, что погрешность двухточечной формулы с возрастанием т при т > 1 уменьшается существенно быстрее, чем погрешность формулы Эйлера-Маклорена.
Из данных, представленных в табл. 4 и из рис. 3, следует, что, например, при т = 7 величина погрешности, полученной по двухточечной формуле интегрирования более, чем на шесть порядков меньше погрешности, полученной по формуле Эйлера-Маклорена для значений п = 1,2,4.
Численный пример 2
Рассмотрим интеграл
I _ [ 2 Лх Л х '
который может быть вычислен аналитически и значение которого равно 1п 2 = 0.69314718....
Для подынтегральной функции в этом случае также существуют производные сколь угодно высокого порядка, которые представляются формулой
(j)
(-1) jj! xj+1 '
(30)
С использованием этой формулы и квадратурных формул, представленных в табл. 1, и формулы (28), получены численные значения интегралов с использованием двухточечного интегрирования и по формуле Эйлера-Маклорена.
Полученные численные данные в виде зависимости погрешности 8(т) при значении параметра п _ 1,2 для обоих методов представлены на рис. 4.
Рис. 4. Сравнение погрешности по формуле двухточечного интегрирования и Эйлера-Маклорена для функции /(х) _ 1/х.
Из графиков, представленных на этом рисунке, видно, что с увеличением т погрешность интегрирования 8т, рассчитанная по двухточечной формуле, монотонно убывает для всех п.
В тоже время погрешность интегрирования, полученная по формуле по Эйлера-Маклорена, с возрастанием т сначала убывает, затем, начиная с некоторого т, начинает резко увеличиваться, т.е. процесс приближения становится расходящимся. Значение критического числа т, начиная с которого погрешность начинает возрастать, зависит от п — числа точек разбиения отрезка интегрирования и монотонно увеличивается с увеличением п.
Такое поведение погрешности формулы по Эйлера-Маклорена для данного интеграла объясняется рассмотрением остаточного члена, модуль которого в соответствии с (26), (27), (30) и с и учетом Ь _ 1, может быть записан в виде
_ 2(2т + 2)! S2m+2 (1 2) |rEm| _ ^--Л2т+2 n2т+3 ' п (1' 2)
(2п П)2
Из того, что первая дробь стремится к бесконечности при т ^ те и ограниченности второй дроби снизу некоторой константой при т ^ те следует, что |гЕт| ^ ^ при т ^ те, т.е. процесс приближения расходится при любом фиксированном п. Однако это не исключает того, что при заданном значении п погрешность может сначала уменьшаться, достигнув некоторого минимума, и только затем начнет неограниченно возрастать, что и показывают результаты численного эксперимента, представленные на рис. 4.
Из формулы (31) следует, что модуля остаточного члена |гЕт| имеет место Лт < |гЕт|, где оценку погрешности снизу Лт можно представить в виде
Ат = Ыф (т),
* (m) = sm& (32)
(2т + 2)! (2п п)2
и М — некоторая константа.
На рис. 5 представлены графики зависимости функции ф(т) для значений параметра п = 1,2,3,4.
Рис. 5. Зависимость функции ф(т) для значений параметра п = 1,2,3,4
Сравнение рис. 4 и рис. 5 показывает, что поведение зависимостей погрешности и ее оценки качественно имеет близкий характер; при этом оценка погрешности, достигая некоторого минимума, также неограниченно увеличивается при т ^ те и минимальные значения т, при которых функция ф(т) принимает минимальные значения, увеличиваются с увеличением п.
Качественный характер поведения функции ф(т) можно исследовать, используя формулу Стирлинга [15, с. 491]
/ п\п .--1 1
"! = Ь)^(1 + УЫ + 288? + -)- (33)
Вводя переменную к = 2т + 2 и ограничиваясь приближенным представлением п! в виде
п! = (=)", 100
которое, однако, позволяет получить аналитическое решение, для функции ф(к) в соответствии с (32) получим следующее представление:
ф (к) _ ' к
2п вп/
Для исследования функции ф(к) дискретного аргумента к используем функцию
ф (х)_ (а )х
непрерывной переменной x и соответствующий аппарат математического анализа.
Дифференцируя функцию ф(x), приравнивая производную нулю в соответствии с достаточным условием экстремума и решая линейное уравнение, получим для точек минимумов x* формулу x* = 2пп. Возвращаясь к переменной m, получим для точек минимума приближения функции ф(m) соотношение
m* = п п — 1.
Из этой формулы следует, что положения минимумов оценки погрешности увеличиваются с увеличением п, что подтверждается результатами расчета, представленными на рис. 5.
Из полученных результатов следует, что формула Эйлера-Маклорена для данного интеграла в случае фиксированного числа точек п приводит к тому, что процесс приближения при m ^ ^ расходится. Однако при этом для фиксированного порядка производной m можно подобрать число точек п так, чтобы приблизить интеграл с заданной точностью. Данный пример подтверждает, что в ряде случаев интегрирования с использованием формулы Эйлера-Маклорена образуется соответствующий ей ряд, который, как отмечено в [10, с. 544], «вообще говоря, расходится».
В тоже время из представленных на рис. 4 результатов видно, что погрешность 8m интегрирования по двухточечной формуле монотонно убывает как с увеличением m, так и с увеличением п. При этом погрешность двухточечного интегрирования так же, как и в примере 1, начиная с m = 2 меньше, чем погрешность результатов, полученных по формуле Эйлера-Маклорена.
Сходимость результатов двухточечного интегрирования для данного примера определяется поведением остаточного члена rlm, модуль которого в соответствии с (21), (30) и с и учетом L = 1, может быть записан в виде
b
Irlm\ = n2m+2^2m+3 > M G (1 2). (34)
Поскольку коэффициенты bm, определенные (14), стремятся к нулю при m ^ <*> (см. например, [16]), то модуль остаточного члена r/m также стремятся к нулю при m ^ <*>, что и обуславливает сходимость результатов двухточечного интегрирования для данной функции.
Заключение
В рамках решения задачи представления определенных интегралов с использованием производных получены составные квадратурные формулы, использующие производные заданного порядка на отрезке интегрирования.
Дано представление для остаточного члена, на основе которого получена оценка погрешности интегрирования, выраженная через производную функции соответствующего порядка.
Рассмотрен пример численного интегрирования функции, для которой известно точное значение интеграла и для которой существуют и ограничены производные сколь угодно высокого порядка. Приведенные численные данные о погрешности интегрировании и ее оценке показывают в этом случае сходимость интегральных приближений к интегралу от заданной функции, как при увеличении числа узловых точек, так и повышении порядка используемых производных.
Дано сравнение метода двухточечного интегрирования с подходом, основанном на использовании формулы Эйлера-Маклорена, которое показало, что для достаточно гладких функций точность двухточечного интегрирования существенно выше, чем по формуле Эйлера-Маклорена. Приведен пример интеграла, для которого его приближения, полученные с использованием формулы Эйлера-Маклорена, расходятся, а по формуле двухточечного интегрирования сходятся и достаточно быстро.
Предложенные составные квадратурные формулы могут использоваться в теоретических исследованиях и практических приложениях для вычисления определенных интегралов.
Конкурирующие интересы. Конфликтов интересов в отношении авторства и публикации нет.
Авторский вклад и ответственность. Автор участвовал в написании статьи и полностью несет ответственность за предоставление окончательной версии статьи в печать.
Список литература/References
[1] Березин И. С., Методы вычислений, Физматлит, М., 1962, 464 с. [Berezin I.S., Metody vychisleniy, Fizmatlit, M., 1962 (in Russian), 464 pp.]
[2] Бахвалов Н. С., Жидков Н.П., Кобельков Г. М., Численные методы, Учеб. пособие, Наука, М., 1987, 600 с. [Bakhvalov N.S., Zhidkov N. P., Kobel'kov G.M., Chislennyye metody, Ucheb. posobiye, Nauka, M., 1987 (in Russian), 600 pp.]
[3] Вержбицкий В.М., Численные методы, ООО Издательский дом Оникс 21 век, М., 2005, 400 с. [Verzhbitskiy V.M., Chislennyye metody, OOO Izdatel'skiy dom Oniks 21 vek, M., 2005 (in Russian), 400 pp.]
[4] Волков Е.А., Численные методы, Учеб. пособие для вузов, Наука, М., 1987, 248 с. [Volkov Ye.A., Chislennyye metody, Ucheb. posobiye dlya vuzov, Nauka, M., 1987 (in Russian), 248 pp.]
[5] Хемминг Р. В., Численные методы для научных работников и инженеров, Наука, М., 1972, 400 с. [Khemming R.V., Chislennyye metody dlya nauchnykh rabotnikov i inzhenerov, Nauka, M., 1972 (in Russian), 400 pp.]
[6] Шахов Ю.Н., Деза Е. И., Численные методы, Учебное пособие, Книжный дом Либ-роком, М., 2010, 248 с. [Shakhov Yu.N., Deza Ye. I., Chislennyye metody, Uchebnoye posobiye, Knizhnyy dom Librokom, M., 2010 (in Russian), 248 pp.]
[7] Калиткин Н.Н., Численные методы, Учебное пособие, БХВ-Петербург, СПб., 2011, 592 с. [Kalitkin N.N., Chislennyye metody, Uchebnoye posobiye, BKHV-Peterburg, SPb., 2011 (in Russian), 592 pp.]
[8] Крылов В. И., Приближенное вычисление интегралов, Наука, М., 1967, 500 с. [Krylov V. I., Priblizhennoye vychisleniye integralov, Nauka, M., 1967 (in Russian), 500 pp.]
[9] Крылов В. И., Квадратурные формулы, Наука, М., 1988, 256 с. [Krylov V. I., Kvadraturnyye formuly, Nauka, M., 1988 (in Russian), 256 pp.]
[10] Фихтенгольц Г. М., Курс дифференциального и интегрального исчисления. Т. 2, Наука, М., 1970, 800 с. [Fikhtengol'ts G. M., Kurs differentsial'nogo i integral'nogo is-chisleniya. V. 2, Nauka, M., 1970 (in Russian), 800 pp.]
[11] Шустов В. В., "О представлении интегралов значениями функции и ее производных на основе использования двухточечных многочленов Эрмита", Теория операторов, комплексный анализ и математическое моделирование, тезисы докладов XIII Межд. научн. конф., Владикавказ, ЮМИ ВНЦ РАН, 2016, 85-87. [Shustov V. V., "O predstavlenii integralov znacheniyami funktsii i yeye proizvodnykh na osnove is-pol'zovaniya dvukhtochechnykh mnogochlenov Ermita", Teoriya operatorov, kompleksnyy analiz i matematicheskoye modelirovaniye, tezisy dokladov XIII Mezhd. nauchn. konf., Vladikavkaz, YUMI VNTS RAN, 2016, 85-87 (in Russian)].
[12] Шустов В. В., "О приближении функций двухточечными интерполяционными многочленами Эрмита", Журнал вычислительной математики и математической физики, 55:7 (2015), 1091-1108. [Shustov V.V., "O priblizhenii funktsiy dvukhtochechnymi interpolyatsionnymi mnogochlenami Ermita", Zhurnal vychislitel'noy matematiki i matem-aticheskoy fiziki, 55:7 (2015), 1091-1108 (in Russian)].
[13] Кудрявцев Л. Д., Математический анализ. Т. 1, Высшая школа, М., 1970 (in Russian), 592 с.
[14] Корн Г., Корн Т., Справочник по математике для научных работников и инженеров, Наука, М., 1984, 832 с. [Korn G., Korn T., Spravochnik po matematike dlya nauchnykh rabotnikov i inzhenerov, Nauka, M., 1984 (in Russian), 832 pp.]
[15] Грехем Р., Кнут Д., Паташник О., Конкретная математика. Основание информатики, Мир, М., 1998, 703 с. [Grekhem R., Knut D., Patashnik O., Konkretnaya matematika. Osnovaniye informatiki, Mir, M., 1998 (in Russian), 703 pp.]
[16] Шустов В. В., "О представлении интегралов значениями функции и ее производных на основе использования двухточечных многочленов Эрмита", Математический форум. (Итоги науки. Юг России). ЮМИ РАН, 11 (2017), 113-122. [Shustov V.V., "O predstavlenii integralov znacheniyami funktsii i yeye proizvodnykh na osnove ispol'zovaniya dvukhtochechnykh mnogochlenov Ermita", Matematicheskiy forum. (Itogi nauki. Yug Rossii). YUMI RAN, 11 (2017), 113-122 (in Russian)].
[17] Шустов В. В., "О представлении определенных интегралов значениями функции и ее производных", Владикавказский математический журнал, 22:2 (2020), 83-97. [Shustov V. V., "O predstavlenii opredelennykh integralov znacheniyami funktsii i yeye proizvodnykh", Vladikavkazskiy matematicheskiy zhurnal, 22:2 (2020), 83-97 (in Russian)].
Vestnik KRAUNC. Fiz.-Mat. Nauki. 2021. vol. 34. no. 1. P. 88-104. TSSN 2079-6641
MSC 41A55 Research Article
On approximation of definite integrals by compound quadrature
formulas using derivatives
V. V. Shustov
State Research Institute of Aviation Systems, 125319, Moscow, Viktorenko st., 7, Russia E-mail: [email protected]
The problem of computing a definite integral of a function for which the values of itself and the set of derivatives up to a given order at the points of the interval of integration are known is considered. Composite quadrature formulas are constructed that use the values of the function and its derivatives up to the m-th order inclusive. A representation of the remainder is obtained, expressed in terms of the derivative of the corresponding order and the number of nodal points. Examples of integration of the given functions with the study of the error and its estimation are given. A comparison is made with the known numerical methods and the Euler-Maclaurin formula, which showed increased accuracy and better convergence of the two-point integration method.
Key words: quadrature formulas using derivatives, compound quadrature formula, remainder of integration, approximation error estimate, Euler-Maclaurin formula
DOT: 10.26117/2079-6641-2021-34-1-88-104
Original article submitted: 20.10.2020 Revision submitted: 03.02.2021
For citation. Shustov V. V. Analysis of the possibility of using Russian geomagnetic stations to calculate the analogue of the Dst-index. Vestnik KRAUNC. Fiz.-mat. nauki. 2021,34: 1,88-104. DOT: 10.26117/2079-6641-2021-34-1-88-104
Competing interests. The author declares that there are no conflicts of interest regarding authorship and publication.
Contribution and Responsibility. The author contributed to this article. The author is solely responsible for providing the final version of the article in print. The final version of the manuscript was approved by the author.
The content is published under the terms of the Creative Commons Attribution 4.0 International License (https://creativecommons.org/licenses/by/4.0/deed.ru)
© Shustov V. V., 2021
Funding. The study was carried out without financial support from foundations.