Алгоритмы численного интегрирования с правилом остановки Текст научной статьи по специальности «Математика»

Известия Тульского государственного университета Естественные науки. 2013. Вып. 3. С. 193-201 = Прикладная математика и информатика =

УДК 511.3

Алгоритмы численного интегрирования с правилом остановки

Н. К. Серегина

Аннотация. В работе рассматриваются алгоритмы численного интегрирования с правилом остановки для параллелепипедальных сеток. Приводятся программы в системе MATHCAD15, реализующие эти алгоритмы.

Ключевые слова: параллелепипедальная сетка, алгоритм с правилом остановки, квадратурные формулы, метод оптимальных коэффициентов, количественная мера качества сетки.

Актуальность и цель работы. Актуальность темы состоит в том, что построение и исследование алгоритмов численного интегрирования с правилом остановки являются важным современным направлением в теории теоретико-числовых методов приближенного анализа (см. [4]).

Цель работы — подробно рассмотреть классический случай теоретикочислового метода Коробова оптимальных коэффициентов интегрирования периодических функций многих переменных, а также дать программную реализацию алгоритмов численного интегрирования с правилом остановки для последовательности концентрических параллелепипедальных сеток.

Как показали Н. Н. Ченцов и И. Ф. Шарыгин вопрос об интегрировании произвольных функций с помощью периодизации задачи сводится к случаю периодических функций. Поэтому рассматриваем только классический случай теоретико-числового метода Коробова — интегрирование периодических функций многих переменных. В работе предлагаются соответствующие программы в системе компьютерной математики МаЛ-еаё15, реализующие данные алгоритмы численного интегрирования.

Параллелепипедальные сетки Н. М. Коробова. В 1959 году профессор Н. М. Коробов предложил новый класс теоретико-числовых сеток

— параллелепипедальные сетки:

Введение

и соответствующие квадратурные формулы с равными весами

/■■■/'*• - NЁ /({£}.....Ш) -

Км [/],

где Км [/] - погрешность квадратурной формулы.

На классе Еа периодических функций с быстро сходящимися рядами Фурье были получены наилучшие результаты

1П«(«-1) N

|Км[/]| << ---тт---- (Н. С. Бахвалов, Н. М. Коробов).

N а

Качество параллелепипедальной сетки. Количественной мерой качества набора коэффициентов а0, а1,..., а3 параллелепипедальной сетки называется величина

нм —^ Ё I! (! - 2{а1г})2- (1)

р к=0 .7=0 4 К р ) /

которая равна приближенному значению интеграла от функции Л,(х) —

— ^ П(1 — 2х7-) по квадратурной формуле с параллелепипедальной

.7=0

сеткой

1 — /1 ... / Вд«х — ^ £ 17 (1 — 2 ■ { }) 2 — ЛрН

■]о -]о р к=0 7=Л I р )/

где Кр[Ь] — погрешность приближенного интегрирования.

Выбор функции Л,(х) и величины Н(р, а) связан с тем, что функция Л,(х)

является граничной функцией класса Еа ^-, Пт) (подробности см. в [10]).

Вычисление количественной меры качества набора коэффициентов. Следующая программа вычисляет по заданному модулю р и набору коэффициентов а — (а0,а1,... ,а3)т значение функции Н (р, а), которая является количественным критерием качества набора коэффициентов для использования в квадратурных формулах с параллелепипедальными сетками:

НТ(р,а) : =

в ^ 1впдЫг(а) — 1

р-1

Е

я ^ ^

р

|_к=о

я

По [1 — 2 •( ^ — Цаат^ ^))]

Например, при р — 101 и а — (1,19, 85)т получим

НТ(р, а) — 1.1030731532962952.

Ясно, что трудоемкость вычисления по этой программе равна О (в ■ р) элементарных арифметических операций.

1

1

т

Вычисление точек сетки. При реализации произведения двух параллелепипедальных сеток естественно иметь программу вычисления точек сетки в виде двумерного массива, что и делает программа БеЬка(р, а):

Setka(p,a) :=

s — length(a) — І, Sp—1,s for k Є 0..p — І forj Є 0..s

0

у-* a a ' k r> і / a a ' k

Sk,j —------------3— — floor \ j-

Эта программа также требует для формирования двумерного массива координат точек сетки О (в ■ р) элементарных арифметических операций. Следующая программа:

Setka1(p,a) :=

s — length(a) — І, Sp—1,s — 0 for k Є І.^ — І forj Є 0..s

r —— Sk—1,j + aj Sk j — if (r < p,r,r — p)

S S '

S — — p

Программа Setka1(p,a) при условии, что выполнены неравенства 1 ^ ^ aj < p (j = O,..., s) обеспечивает более точное вычисление координат точек при больших значениях p, так как операция деления выполняется один раз над всем массивом, а при вычислении числителей координат все вычисления ведутся с числами в диапазоне от О до 2p и не происходит потери точности из-за вычитания двух плавающих чисел из разных диапазонов значений, которое встречается в первой программе.

Вычисление количественной меры качества сетки. При наличии массива координат точек параллелепипедальной сетки естественно иметь программу, которая вычисляет значение функции H непосредственно по координатам точек. Таким образом получаем программу HS(S):

HS(S) :=

p — rows(S), s — cols(S) — І -p— 1 -E

k=o

R—

R

3s+1

p

П (І — 2 ' Sk,j)

j=o

Отметим, что, несмотря на то, что МАТНСАБ15 является интерпретатором, перечисленные выше программы работают достаточно эффективно и требуют не более одной минуты времени на персональном компьютере средней производительности при р — 5000000.

S

2-

1. Оптимальные параллелепипедальные сетки

1.1. Оптимальные коэффициенты сетки. Покоординатный алгоритм. Переходя к вопросу о вычислении оптимальных коэффициентов по модулю р, прежде всего рассмотрим покоординатный алгоритм вычисления, который требует для своей реализации О (в2 ■ р2) элементарных арифметических операций:

НТ 1(р,в) :=

а5-х —— 1, Ьо —— 1, ао —— 1 $огз € 1..в — 1 Я — 3,+1

for с € 1..р — 1

Ь, — с

3^ + 1

р-1

Е

к=0

==о ^ — 2( >_£ — ,юог( ь^))]

с, Я ■

гЦ г < Я

Ь, — а,

1.2. Обоснование покоординатного алгоритма. Необходимо отметить, что первоначально Н. М. Коробов свой алгоритм обосновал только для простых модулей. Позднее он его модифицировал для составных модулей равных произведению двух простых, что позволило сократить трудоемкость алгоритма вычисления оптимальных коэффициентов с О (Ж2) до О (Ж). Хотя в литературе не встречается обоснование этого алгоритма для произвольного составного модуля, но не трудно понять, что такое обоснование, но без оценок качества, нетрудно дать.

Действительно, из результатов разных авторов (например, см. [3]) следует существование оптимальных коэффициентов по любому составному модулю. Кроме того, известно, что в-мерный набор оптимальных коэффициентов всегда можно дополнить до (в + 1)-мерного набора. Поэтому, если среди всех (в + 1)-мерных наборов с фиксированными первыми в коэффициентами взять набор с минимальным значением функции Н — количественной меры качества оптимальных коэффициентов, то этот набор и будет оптимальным.

1.3. Произведение сеток. Пусть х = (ж1,...,ж8) € М5, тогда {х} = = ({х1},..., {х5}) — дробная часть х.

Для любых двух сеток М1 и М2 их произведение есть сетка

М = М1 ■ М2 = {{х + у}|ж € М1,у € М2}.

Различные аспекты теории произведения обобщенных параллепипедальных сеток рассматриваются в работе [5].

При построении алгоритма Л. П. Добровольской поиска оптимальных коэффициентов существенно используется представление параллелепипе-

2

Г ——

Р

,

а

дальной сетки по составному модулю как произведение сеток по взаимно простым модулям.

1.4. Модифицированный алгоритм Л. П. Добровольской.

Следующая программа по заданной сетке 5 из Р точек, представленной двумерным массивом, находит оптимальную параллелепипедальную сетку из р точек, такую, что их произведение будет оптимальной параллелепипедальной сеткой из Р ■ р точек. Важной особенностью данной программной реализацией является тот факт, что массив с результирующей сеткой представляет из себя р подмассивов, каждый из которых суть модифицированная исходная параллелепипедальная сетка. Эта особенность потом будет использована при реализации алгоритма численного интегрирования с правилом остановки.

Данный алгоритм является модификацией алгоритма Л. П. Добровольской.

Здесь вместо очень сложной системы целевых функций используется количественная мера качества сетки Н(5):

HTS(S,p, s) :=

as—1 — І, bo — І, ao — І, P — rows(S) for j Є І.. s - І R — 3j+1 for с Є І..р — І

3j+1 p • p

P—1 p—1

E E

.n = 0 k = 0

— floor (^Sn,i + ^pr))]

aj — c, R — r if r < R

SP • p—1,s —1 — 0

for k Є 1..p — 1 for n Є 0.. P — І for j Є 0..s — 1

Sp • k+n,j — Snj + apk — floor (Snj +

bj —c

r —

2

S

Программа НТ5(5, р, 8) требует для своего выполнения память объемом О (8 ■ Р ■ р) и О (з2 ■ Р ■ р2) элементарных арифметических операций.

1.5. Об алгоритме Л. П. Добровольской. Алгоритм Л. П. Добровольской основан на двух идеях.

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

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

1.6. Допустимая последовательность простых. Последовательность ро, ..., ри называется допустимой последовательностью простых, если для Р = ро ■ ... ■ ри выполнены соотношения

Р

ри ■ ри-1 ■ ... ■ р,'+1 ■ р2 < 2Т-2 (1 ^ ^ к)-

Для этого достаточно, чтобы

рз <

р^-1 ■ ... ■ р1 ■ ро 2з—2

Если

рз—1 ■ ... ■ р1 ■ ро ^ ^ рз-1 ■ ... ■ р1 ■ ро

-—2з—--------< рз <—2—2---------- (1 < 3 <

то получаем монотонную последовательность простых ро < ... < ри.

1.7. Программная реализация модифицированного алгоритма Л. П. Добровольской. Следующая программа НТ51(рр, в) по заданному вектору рр = (ро, ...,ри)Т, где ро, ..., ри допустимая последовательность простых, строит оптимальную параллелепипедальную сетку из Р = ро ■ ... ■ ри точек:

НТЯ1(рр, в) := к — 1впдИг(рр) — 1

а — НТ 1(ррк, в), Я — ЯеЬка(ррк,а) forj € к — 1,к — 2..0 Я — НТЯ(Я,рр, ,в)

Я

1.8. Трудоемкость алгоритма. Нетрудно видеть, что количество операций требуемых для выполнения этой программы для допустимой последовательности простых есть

О (в2(р| + ри ■ р|-1 + ... + ри ■ ри-1 ■ ... ■ р1 ■ ро)) =

= О в

(в2 (Е 2—2 + Р ■ ро)) = О (в2 ■ Р ■ (4+ ро)) .

2. Алгоритм интегрирования с правилом остановки

2.1. Тестовая функция. В качестве тестирующей предлагается функция

__ 0^+1

П(1 - 2 ■ Хз )

3=о

которая реализуется программой

f (х) :=

в — 1впдШ(х) — 1 в

38+1 • П (1 — 2 • х,)

,=о

/ ... / (х)^жо...^ж5 = 1,

оо

2

2.2. Мультипликативная дискретная дисперсия по произведению параллелепипедальных сеток. Пусть Мі и М2 — две параллелепипедаль-ные сетки и М = Мі ■ М2 — их произведение.

Мультипликативной дискретной дисперсией по произведению параллелепипедальных сеток М для периодической функции / (ж) называется величина

2.3. Численное интегрирование с правилом остановки.

Следующая ниже программа реализует алгоритм численного интегрирования с правилом остановки. Предполагается, что сетка 5 вычислена по алгоритму Добровольской. Тогда пользуясь тем, что двумерный массив, задающий сетку, специально организован, данная программа вычисляет приближенное значение интеграла от функции /(х), наращивая количество точек от первой сетки из ри точек, до последней сетки, если это будет необходимо, из Р точек.

В качестве правила остановки берется малость величины мультипликативной дискретной дисперсии:

IntPo(pp, S, f, є) :=

s — cols(S) — 1, k — length(pp) — 1

x — submatrix(S, О, О, О, s)T, М — f (x),m — М, Int — О for n Є О..ppk — І

x — submatrix(S, n, n, О, s)T, r — f(x)

М — r if r > М m — r if r < m otherwise Int — Int + r Int — ,D — 2 • є,j — k — 1, P — ppk

pp k k

while (D > є) • (j ^ О)

D — Int2 for l Є І..ppj — І Ir — О

for n Є О-.ppk — 1

Pn — P • l + n,x — submatrix(S, Pn, Pn, О, s)T

r — f (x)

М — r if r > М

m — r if r < m otherwise

Ir — Ir + r

Ir

p ,D — D + Ir2, Int — Int + Ir

Int — ^ , D —

--j ’

(Int D m М)

— Int2, P — P • ppj ,j — j — 1

ppj

Программа в качестве результата выдает приближенное интеграла, величину мультипликативной дискретной дисперсии, минимальное и максимальное значения функции на сетке.

значение а также

Работа может быть использована в рамках учебного процесса при изучении квадратурных формул в курсах вычислительной математики и численных методов, а также в курсе по выбору, посвященном изучению теоретико-числовых методов в приближенном анализе

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

1. Бабенко К.И. Основы численного анализа. М.: Наука, 1986.

2. Бахвалов Н.С. О приближенном вычислении кратных интегралов // Вестник МГУ. Сер.1. Математика. 1959. №4. С.3-18.

3. Бочарова (Добровольская) Л.П. Алгоритмы поиска оптимальных коэффициентов // Чебышевский сборник. Тула: Изд-во ТГПУ им. Л. Н. Толстого. 2007. Т. 8. Вып. 1(21). С. 4-109.

4. Добровольская Л.П., Добровольский Н.М., Симонов А.С. О погрешности приближенного интегрирования по модифицированным сеткам // Чебышевский сборник. Тула: Изд-во ТГПУ им. Л.Н. Толстого. 2008. Т. 9. Вып. 1(25). С. 185-223.

5. Добровольский М.Н., Добровольский Н.М., Киселева О.В. О произведении обобщенных параллелепипедальных сеток целочисленных решёток // Чебышевский сборник Тула: Изд-во ТГПУ им. Л.Н. Толстого. 2002. Т. 3. Вып. 2(4). С. 43-59.

6. Коробов Н.М. Вычисление кратных интегралов методом оптимальных коэффициентов // Вестник МГУ. Сер. 1. Математика. 1959. № 4. С. 19-25.

7. Коробов Н.М. О приближенном вычислении кратных интегралов // ДАН СССР. 1959. Т. 124. № 6. С.1207-1210.

8. Коробов Н.М. Свойства и вычисление оптимальных коэффициентов // ДАН СССР. 1960. Т. 132. № 5. С.1009-1012.

9. Коробов Н.М. Теоретико-числовые методы в приближенном анализе. М.: Физматгиз, 1963.

10. Коробов Н.М. Теоретико-числовые методы в приближенном анализе. М.: МЦНМО, 2004.

11. Локуциевский О.В., Гавриков М.Б. Начала численного анализа. М.: ТОО Янус, 1995.

12. Ребров Е.Д. Алгоритм Добровольской и численное интегрирование с правилом остановки // Чебышевский сборник. Тула: Изд-во ТГПУ им. Л.Н. Толстого. 2009. Т. 10. Вып. 1(29). С. 65-77.

13. Огородничук Н.К, Ребров Е.Д. Об алгоритме численного интегрирования с правилом остановки // Алгебра и теория чисел: современные проблемы и приложения: матер. 7 Международной конф. Тула: Изд-во ТГПУ им. Л.Н. Толстого, 2010. С. 153-158.

14. Огородничук Н.К., Ребров Е.Д. ПОИВС ТМК: Алгоритмы интегрирования с правилом остановки // Многомасштабное моделирование структур и нанотехнологии матер. Международной научно-практической конф., посвященной 190-летию со дня рождения академика Пафнутия Львовича Чебышёва, 100-летию со дня рождения академика Сергея Васильевича Вонсовского и 80-летию со дня рождения члена-корреспондента Виктора

Анатольевича Буравихина, Тула: Изд-во ТГПУ им. Л.Н. Толстого, 2011. С. 153-158.

15. Algorithms for computing optimal coefficients / N.M. Dobrovolskiy [et al.] // Computer Algebra and Information Technology: book of abstracts of the International scientific conf. Odessa, August 20-26. 2012. P. 22-24.

16. Некоторые вопросы теоретико-числового метода в приближенном анализе / Л.П. Добровольская [и др.] // Алгебра и теория чисел: современные проблемы и приложения: тр. X Международной конф. Ученые записки Орловского государственного университета. 2012. № 6. Ч. 2. С. 90-98.

Серегина Надежда Константиновна (nadj_nadj@mail.ru), аспирант, кафедра алгебры, математического анализа и геометрии, Тульский государственный педагогический университет имени Л.Н. Толстого, Тула.

Algorithms for numerical integration with stopping rule

N. K. Seregina

Abstract. The paper deals with algorithms for numerical integration of the stopping rule for Parallelepipedal nets. We present a program in MATHCAD15 that implement these algorithms.

Keywords: parallelepipedal net, algorithm with a stopping rule, quadrature formula, method of optimal coefficients, quantitative measure of the quality of the grid.

Seregina Nadegda (nadj_nadj@mail.ru), postgraduate ctudent, department of mathematical analysis, algebra and geometry, Leo Tolstoy Tula State Pedagogical University.

Поступила 23.08.2013