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

УДК 519.6+517.9

КВАДРАТУРНЫЕ ФОРМУЛЫ ВЫСОКОГО ПОРЯДКА АППРОКСИМАЦИИ

Д-А. Силаев

В работе предлагается метод построения квадратурной формулы высокого порядка аппроксимации для широкого класса областей, основанный на приближении гладкой функции на плоскости полулокальным сглаживающим сплайном или Б-сплайном. Полулокальные сглаживающие сплайны были введены Д.А. Силаевым. Ранее рассматривались и применялись сплайны 3-й и 5-й степени. Настоящая работа посвящена использованию Б-сплайнов более высоких степеней. Появление устойчивых Б-сплайнов класса С0 (только непрерывных), состоящих из полиномов высокой степени п (п = 9,10) позволило получить квадратурные формулы 10-го и 11-го порядков аппроксимации. Предполагается, что интегрируемая функция принадлежит классу Ср (р = 10,11) в несколько большей области, чем исходная область, по которой ведется интегрирование. Предполагается также, что граница области задана параметрически, что позволяет с высокой степенью точности учесть границу области. Подобный подход возможен и для построения кубатурных формул.

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

Введение

Теория квадратурных формул направлена на получение приближенных формул для вычисления интеграла, максимально точных при наименьшем числе узлов [1-4]. Составные квадратурные формулы (формула трапеций, Симпсона и т.п.) можно интерпретировать как формулы, полученные с помощью приближения интегрируемой функции интерполяционным сплайном типа Лагранжа (локальными сплайнами). Использование глобальных сплайнов также приводит к квадратурным формулам [5, 6]. Функции многих переменных могут быть приближены суммой произведений функций одной переменной [7], что приводит к построению квадратурных и кубатурных формул. Однако, существенным недостатком использования глобальных сплайнов является отсутствие удобных алгоритмов их построения для случая высоких степеней (выше 3, 5). Остается большой трудностью аппроксимация интегрируемой функции в окрестности границы области [8, 9].

На плоскости рассматривается ограниченная область П с границ ей 7 = д П. Предполагается, что граница задана параметрически

7 : [х(г),у(г), \г Є [а,в]}, где х,у Є

С!+е

- заданные периодические функции, т.е. х(а) = х(в), у(а) = у(в) , первые производные функций по х, у удовлетворяют условию Гельдера с порядком є > 0. В несколько большей области рассматривается гладкая функция / Є Сга+1, т.е. она имеет непрерывные ограниченные п + 1 частные производные. Квадратурная формула имеет вид:

JJ /(х,у)йхйу\ = си/(Р^ + 0{Нп+1), (1)

где Ск - теса, Рк - узлы квадратурной формулы, Н - шаг разбиения. Здесь к - номер полинома, составляющего сплайн, отвечающего точке Рк, а Ск - интеграл по области О от указанного полинома [10, 11]. Дополнительно будем предполагать, что исходная функция такова, что она определена в несколько большей области О$, содержащую область О с сохранением класса и нормы. Для простоты можно считать, что таковой областью является круг. Проблема, с которой мы здесь сталкиваемся, состоит из двух моментов. Во-первых, как унифицировать вычисление большого числа таких интегралов по заданной области. Во-вторых,

О

О

границе области 7 = дО.

1. Одномерный 5-сплайн класса Су

Рассмотрим на отрезке [а, Ь] равномерную сетку {хк}к=К ^ хк = а + кН, Н = (Ь — а)/К — шаг сетки. Разобьём отрезок [а, Ь] на группы, доя этого введём на [а, Ь] ещё одну равномерную сетку {^тУ^о, £1™ = а + Ш7 Н = тН, т € N . Здесь т фиксировано, индекс т в обозначении в дальнейшем будем опускать. Таким образом, переходя из одной группы в другую, мы осуществляем сдвиг системы координат и рассматриваем каждый 1-й полином на отрезке [0, Н] . Пусть значения приближаемой функции на этой сетке у = (уо,у\,..., ук) € И,к+1. Обозначим через

П |

и : и(х) = а0 + а1х ... + архр + ^ aixг

г=р+1 J

множество полиномов степени п с фиксированными коэффициентами ао,а1,... ,ар. Рассмотрим функционал

м

Ф (и) =^2(и(& + кН) — Ут1+к)2 . к=0

В классе РПП ищется такой пол ином д1, который минимизирует функционал

м

ф1(и) = ^2(и(& + кН) — ут1+к )2 —► шт(ар+1, ...,ап) к=0

и удовлетворяет следующим условиям:

ао = д-1(Н),а[ = д—1(Н),...,аР = Уд^1(Н) при I = 1,...,ь — 1. (2)

у.

В случае периодического 5—сплайна здесь при I = 0 выполнено д-1(Н) = до-1 (Н). Так

как

а0 = д1 (0),а1 = д'М, ...,а1т = тд1т)(0), при т = 0,1,... ,р,

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

(р)

уо,у0,..., У0Г1- Можно предполагать, что значения заданной функции ук известны с некоторой точностью, например, они есть результаты каких-либо измерений. Будем предпола-

Н

1В случае если функция задана таблицей, то у0, у'0,..., удР\ можно вычислить с помощью формул численного дифференцирования высокого порядка аппроксимации (см. [12]).

будем предполагать, что если функция f € Сп+1[а,Ъ] задана в узлах равномерной сетки хи = а + кН, к = 0,1,..., К своими значениями у&, то \уи — f (хи)| < ОН(п+1'). Здесь Ь - число групп, на которые разбита исходная таблица значений приближаемой функции, или число полиномов, составляющих сплайн. Кроме того, здесь М + 1 - количество точек осреднения, т + 1 - количество точек, входящих в область определения 1-то полинома д^ £I - точка привязки полинома д^ М — т + 1 - число таких точек, значения которых участвуют при определении двух соседних полиномов, составляющих 5-сплайн, М > т + 1 ([12]).

Определение 1. Б-сплайном назовем функцию Бт,м(х), которая совпадает с полиномом дI (х) на каждом отрезке £I < х < £1+1.

Система линейных алгебраических уравнений, которой должны удовлетворять коэффициенты полиномов Б-сплайна, состоит из уравнений двух видов: а) уравнений склейки для каждой пары последовательных полиномов (2); б) уравнений для определения коэффициентов при старших степенях полиномов по коэффициентам при младших степенях. Сделаем замену переменных с = аН, г = 0,1,... ,п.

Обозначим через

м

м

Б, = £ кз, Р) = £ ут1+ккр+, і = 1,...,и — р ,СП =

и

р\(и — р)!

(3)

к=0 к=0

Здесь I = 0, . . . , Ь — 1 - номер полинома, причем если I = 0, то в периодическом случае выражение а1— 1 означает а^—1. В дальнейшем, если это не вызовет путаницы, мы будем опускать волну над переменными а\. Запишем систему для определения коэффициентов полиномов в матричной форме. Для этого обозначим через

БР+1 . . Б2р+1 Б2р+2 . . Бп+р+1

А0 = Бп ■ Бп+р , А = Бп+р+1 . . Б2п

1 т 2 т2 . . . тр тр+1 тр+2 . . . тп

0 О 0 1 2т.. . ртр-1 , В1 = (р + 1)тр (р + 2)тр+1 . . . итп-1

0 0 0 .. .1 СР+1т 2 .. 2 . СРтп—Р

Здесь прямоугольные матрицы Ао и В1 имеют размерности (п — р) х (р + 1) и (р + 1) х (п — р) соответственно, размерности квадратных матриц А1 и Во (п — р) х (п — р) и (р + 1) х (р + 1). Пусть, кроме того,

I Р1 \

Р1 Р2

*0 =

V р1 /

\ п—р /

( а0 \

\аР

( аР

XI =

+1

\

ар+2

\ аП /

ГД6

I = 0,1,...,Ь — 1. (4)

Тогда уравнения а) склейки для каждой пары последовательных полиномов (2) примут вид:

ВоХО + ВХ1 = Х0+1 , (5)

а уравнения б) для определения коэффициентов при старших степенях полиномов по коэффициентам при младших степенях!

лх0 + Аіхі = р

(6)

а

1

2. Существование и единственность 5-сплайна класса Ср

Предположим, что т и М таковы, что матрица А имеет обратную. Тогда из (6) получаем, что

Х[ = А—1Р1 — АХО , (7)

где А = А—1 Ао. Подставим выражение для X1 в (5). Тогда получим рекуррентное соотношение, связывающее р + 1 младших коэффициентов I + 1 полинома через р + 1 младших коэффициентов I полинома:

Х0+1 = иХ10 + Ф1 , (8)

где Ф1 = В1 А—1Р1, матрица устойчивости и = Во — В1А— 1Ао имеет размерность (р + 1) х (р + 1)

Рассмотрим сначала непериодический случай. Зададим начальный вектор

где значения производных, ВХОДЯЩИХ В х0 , могут быть вычислены приближенно с высокой степенью точности с помощью формул численного дифференцирования. Пользуясь формулами (7),(8), последовательно находим Х°,Х0,... ,Х^—1. Тем самым все коэффициенты полиномов, составляющих сплайн, однозначно определены.

Теорема 1. Пусть числа т, М,р, и таковы, что (1еЛА1 = 0 . Тогда для любой функции f (х), заданной на отрезке [а, Ъ\ своими значениями ук в точках Хк = а + кН, Н = (Ъ — а)/К и начального вектора Х[° существует единственный непериодический сплайн Б'^м[у\(х) р

В периодическом случае применяя рекуррентную формулу (8) Ь — 1 раз, получим:

ь

Х0 = ХЬ = ихЬ—1 + ФЬ—1 = и (ихЬ—2 + фь—2) + фь—1 = ... = и ьх0 + ^ и Ь—*ф*—1,

в=1

откуда

х0 = (е — и ь)—1 ^ и ь—афа—1

в=1

Затем последовательно находим Хо,Хо,... ,Х^ 1. Тем самым все коэффициенты полиномов, составляющих периодический сплайн, однозначно определены.

Теорема 2. Пусть числа т, М,р,п таковы, что (1еЛА1 = 0 и собственные числа матрицы и не равны корню степени Ь из единицы (здесь Ь - число полиномов, составляющих сплайн). Тогда, для, любой функции f (х), заданной на отрезке [а, Ь] своими значениями уи в точках хи = а + кН, Н = (Ь — а)/К существует единственный периодический сплайн Б'Пт м[у](х) класса (У.

3. Устойчивость и сходимость 5-сплайна кла сса Ср

Теорема 3. Пусть периодическая функция f (х) € Сп+1 [а,Ь], и пусть выполнено условие

\1'(хи) — Ук \ < ОоНп+1+£, е > 0 . (9)

Пусть, кроме того, числа т, М,р,п таковы, что (1еЛА1 = 0, и собственные значения мати

М < 1, г = 1,...,р + 1. (10)

Тогда, периодический сплайн БП м € Ср[а, Ь] с узлами на, равномерной сетке имеет, дефект п — р, и для х € [а, Ь] справедливы следующие оценки:

Ц-(г)(х) — ^БП,м(х)\< СгНп+1~Т, для, г = 0,1,...,п; (11)

х = & при г = р + 1,... ,п; в этом случае ф(г)(&I) = ф(г)(&I + 0) .

Аналогичные утверждения справедливы и для непериодического случая (см. [12]).

и

йег(и — ХЕ) = 0. (12)

Для случая малых значений М (при 3 < М < 20) в результате расчета были получены значения собственных чисел матрицы и. Как показано в случаях п = 3 и п = 5, для обеспечения условия устойчивости, т.е. выполнения неравенства (10), необходимо перекрывание. Это означает, что имеются такие элементы исходной таблицы значений функции, которые участвуют в определении коэффициентов не менее двух соседних полиномов, составляющих сплайн. Если перекрывание достаточно большое, то это в ряде случаев является и достаточным условием [13]. На практике наиболее употребительными являются те сплайны, для

М

тМ

что в случае р = 0 и п = М, 1 < т < М — 1 матрица и есть число, равное 0. В этих случаях Б—

построения формул численного интегрирования (квадратурных, кубатурных).

4. Фундаментальный 5-сплайн

Фундаментальный периодический Б-сплайн В^ (х) - это Б-сплайн, построенный по данным у = (уо, У1,..., ук) € Кк+1 вида: {у1 = 5^\г,] =0,1,..., К} 5г] - символ Кронекера.

Легко видеть, что линейная комбинация

к

Б(х) = Т1 уз В (х)

з=о

является Б-сплайном, приближающим данные {уг\г = 0,1,...,К}. Непериодические фундаментальные сплайны дополняются сплайнами с начальными условиями у'о,у^,..., у(р)(0), принимающими значения 0 или 1.

5. Двумерный полулокальный сглаживающий сплайн класса Ср

5.1. Построение ф — г — 5-сплайна на круге

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

{фг = гН1\г = 0,1,..., К1}, {Фи = кН1\к = 0,1,... ,Ь1},Н1 = т1Н1, К1 = т1Ь1,К1Н1 = 2п,

{гз = НЦ = 0,1,..., К2}, {Е[ = 1Н\1 = 0,1,..., Ь2}, Н2 = т2Н2, К2 = т2Ь2, К2Н2 = 1.

Будем строить аппроксимацию функции f (ф, г) на круге при условии, что функция f имеет п+1 производных по переменным г и ф, то есть f € С(п+1)[0,1] х [0,2п\. Пусть {угз = f (фг, Г]), г = 0,1,..., К1, ] =0,1,..., К2} - значения в узлах сетки, по которым будет

проводиться аппроксимация. При каждом ] = 0,1,..., К2 построим периодический Б-сплайн Б](ф) на отрезке [0, 2п] по заданным {уц\г = 0,1,..., К\}. Каждый из этих сплайнов аппроксимирует функцию f (ф,г]) на окружности с радиусом г], причём в силу теоремы 3 о сходимости

< CH1+1 ^, ц = 0,,ф € [0, 2п].

йф^^7 j

Далее фиксируем произвольное ф € [0, 2п]. Рассмотрим набор {zj = Sj(ф)\j = 1,..., K2, zo = У00} У00 ~ значение функции f в начале координат. Также обозначим через z0,..., z0p) значения, получаемые по некоторому алгоритму по набору {zj}, которые приближают fr(ф, r)\r=o, • • •, frPP(ф, r)\r=o с порядками n,..., (n + 1 — p) соответственно (например, с помощью формул численного дифференцирования высокого порядка, см. [12]). По набору {zj} и z'0,..., z0p) строим Sp(r) - непериодический S-сплайн на отрезке [0,1]. Будем считать, что m2 < M2Z* ■ Это гарантирует, что собственные значения матрицы U по модулю меньше единицы. Тогда, построенный для ф сплайн S^(r) будет аппроксимировать функцию f (ф, г) при r € [0,1].

— r S( , r)

бых ф и r определяется по следующему алгоритм у: по набору {zj = Sj (ф)\ j = 1,... ,K2, z0 = y00}, z0,..., z0P строим Sv(r), затем полагаем S( ф,r) = Sv(r), другими словами

S (Ф ,r) = {Sv(r)\{zj = Sj(Ф )\ j = 1,...,K2, zo = yoo}}. rn

принадлежащей сетке, то есть r = Rj. При r = Rj определим производную следующим образом:

д v д v

drvS(ф,r) = drvS(ф,r + 0), ^ = 0 1,...,p.

Определение 3. Назовём производной порядка ц по ф от, ф — r-сплайна (1 = 1,... ,n) функцию S (ф ,r) на единичном круге, которая равна ф — r-сплайну, построенному по набору

d^ йф

{zj = ^7Л Sj (ф ^ j = 1,...,K2,zo = У00}.

Как и в случае с производной по г, под производной по ф в точках ф = Фи понимается значение в точке ф = Фи + 0. Наконец, можно ввести понятие смешанной производной.

Определение 4. Под смешанной производной Б ( ф ,г) понимается производная по-

рядка, ц по г от производной порядка V по ф от, Б(ф,г), где производные трактуются согласно определениям 2 и 3.

5.2. Сходимость ф — г — 5-сплайна

Обозначим Н = тах(Н1, Н2).

Теорема 4. Пусть т1/М1 < (*, т2/М2 <(* и f € С(п+1)[0,2п] х [0,1]. Тогда, для,

— г Б( , г)

^(ф,т) — /(ф,т)

< С^Пп+1-^-1У, где /I, V > 0, 0 < ц + V < п. (13)

дт^дфУ дт^дф'

Аналогично можно ввести понятие и т — ф-силайна [14].

5.3. Получение б'-сплайна на круге как явной функции двух переменных

Будем обозначать фундаментальные сплайны по ф как Сг(ф), а фундаментальные сплайны по аргументу т как Dj (т).

Б(ф,т) = {Sф(т)\{Zj = Sj(ф )| з = 1,...,К2,го = уоо}} = 8ф(т).

В свою очередь,

К К2 Кг-1 Кг-1 К2

Бф(т) = ЁzjDj(т) = Dj(т) ^ УVСг(ф) = Ё 5^'С^ф)Dj(т) (14)

j=0 j=0 г=0 г=0 j=0

(здесь ограничиваемся случаем р = 0; в общем случае добавляются слагаемые z,jDjo(т) +

... + ZjpЛ|Djp(т), где фундаментальный сплайн Djp отвечает нулевому набору Zj и х(Р = 1). Предпоследнее равенство следует из определения набора {zj = Sj(ф)} и разложения по фундаментальным сплайнам

К-г-1

Sj(ф) = ^2 у^С^ф)-

г=0

Теперь рассмотрим укрупненную сетку на окружности {Ф* = кИ1\ к = 0,1,... ,Ь\}, где Н1 = т^1 и {К, = 1И2\ I = 0,1,... ,£2}, где И = т2Ь,2■ Рассмотрим вид S-cплaйнa в некотором произвольном секторе этой сетки: ф = кН1 + ф, т = 1Н2 + г, где \ф\ < И1 и \т\ < Н2- В этом секторе фундаментальные S-cплaйны согласно определению представляются в п

пп

С>(ф)=^1 4* Ф, Dj (т) = £ 4^, ф.

р=0 д=0

Подставляя эти выражения в формулу (14) для функции S (ф, т) и меняя порядок суммиро-вэния, получим!

Кг-1 К2 п п п п Кг-1 К2 п п

S(Ф, т)=Т, Е у«И 4*г Е = Е Е фф(Е Е y^j4*4) = Е Е <г"*1■

г=0 j=0 р=0 д=0 р=0 д=0 г=0 j=0 р=0 д=0

Таким образом, показано, что на каждом произвольном секторе функция 5У(ф,т) представ-п

п п К\-1 К2

S(ф,т) = ^фРф:'1, где аР1<1 = £ YlУijсРк^ (15)

р=0д=0 г=0 j=0

или сплайн-функцию двух переменных [14]. Заметим, что в выражения для коэффицие-тов ар1д входят значения всех усодержащихся в круге. Аналогичные выражения можно получить для всех многомерных областей, представляющих собой тензорные произведения одномерных, например, для прямоугольника и тора.

Представление сплайна на круге в виде разложения по одномерным фундаментальным сплайнам (14) позволяет определить понятие смешанной производной для двумерного сплайна.

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

Эту формулу можно рассматривать как формулу численного дифференцирования, основанную на приближении двумерной функции полулокальным сглаживающим сплайном. Всё то же самое верно и для г — ф-сплайна.

5.4. Получение квадратурных формул для одномерных интегралов

Подставим выражение Б-сплайна через фундаментальные сплайны в интеграл:

- искомые коэффициенты квадратуры. Здесь а1?™ - з-й коэффициент т-го полинома в г-м фундаментальном сплайне (т.е. построенном по набору данных {у^ = 5^| к = 0,1,..., К}, где 5^- символ Кронекера). Заметим, что в непериодическом случае указанные фундамен-

/ (р)

тальные сплайны дополняются сплаинами с начальными условиями у0,... ,у0 , принимающими значения 0 или 1. Эти формулы имеют (п + 1)-й порядок аппроксимации.

5.5. Получение квадратурных формул для двумерных интегралов на круге К

Подставим в интеграл по единичному кругу выражение ф — г-сплайна в виде (14):

Определение 5. Под смешанной производной двумерного сплайна

Б(ф, г), где 0 <

^ + V < п, понимается следующая конечная сумма,

к -1 к

(16)

ГД6

Кг-1 К2

Б(ф,г)(Ш= / / Б(ф,г)г(г(ф = ^ ''^Угз С(ф)йф / Б](г)г(г.

г=0 ]=0

(17)

,1 ^2+1 ,.£д + 1 ^2 + 1 ,.И2 П

(1^ = б/(т)тйт = У'' / тБ/(т)йт = У'' / (и + £3) У^ Ъ^п1 йи =

«=о «=0 "'0 1=0

^2+1 . н2 п ь2+1 п , 1 е .

= £ Е Ъ13(пд+1 + Нп1 )йп = £ £ ^Я1+2 ( — + — ) .

«=0 ■J0 1=0 «=0 1=0 'Ч + 4 + '

Здесь а™ш Ъд3 - к-й и д-й коэффициенты ш-го и з-го полиномов в г-м и И-м фундаментальном периодическом сплайне на окружности [0, 2п] и неперпдпческом сплайне на отрезке [0,1 ] соответственно. Здесь Н1 = 2п/Ь1, Н2 = 1 /Ь2■ Фундаментальный периодический сплайн Сг(ф) строится по набору данных {ук = Ьгк| к = 0,1,..., К\}. Непериодический фундаментальный сплайн б/ (т) строится по набору данных {ук = Ь^к | к = 0,1,..., К2; у0,..., у0^}, гДе 5гк~ символ Кронекера, у0. ,у0Р принимают значения либо 0, либо 1.

5.6. Квадратурные формулы для двумерных односвязных областей

На плоскости рассматривается ограниченная область П с границ ей 7 = дП, где 7 -замкнутая самонепересекающаяся кусочно-гладкая кривая. Предполагается, что граница задана параметрически: {^ = {х(Ь) ,у(Ь)}\Ь € [а, в]^, где х,у € С1+£ - заданные периодические функции, т.е. х(а) = х(в),у(а) = у(в), первые производные функций х,у удовлетворяют условию Гельдера с порядком е > 0 (быть может за исключением отдельных точек). Будем предполагать также, что функция / определен а и (п + 1) раз непрерывно дифференцируема в несколько большей области П$. Для простоты будем считать, что область П есть круг К радиуса К.

Построение аппроксимирующей сетки будем производить следующим образом. Поместим область П в круг К радиуса К и введём полярную систему координат, связанную с центром круга. Будем рассматривать в круге радиуса К полярные сетки:

{фг = гН1\г = 0,1,... ,К1}, {Фк = кН1 ,к = 0,1,... ,Ь1},

{тj = НИ = 0,1,..., К2}, {К = 1Н2,1 = 0,1,..., Ь2}, (18)

Н1 = ш1Н1, К1 = ш1Ь1,К1Н1 = 2п, Н2 = ш2Н2, К2 = ш2Ь2, К2Н2 = К.

Пусть пгj = /(фг, Tj) - сужение функции / на равномерную сетку (18). По таблице значений иг;; строим полулокадьный сглаживаю щий сплайн Б (ф, т), состоящий из поли номов п—ой степени, например, т — ф-сплайн, определенный на всем круге К. Из оценки (13) следует, что Б аппроксимирует функцию / с порядком 0(Н(а+1'))1 где Н = шах(Н1,Н2) в области П$. Подставим в интеграл по области П выражение для т — ф-сплайна в виде:

гг К1—1 К2

Б(ф,т)йП = Б(ф,т)тйтйф = У'' У^с4'у^, (19)

п У./п г=0 j=о

где

сг = // Ст,(ф)Б/(т)тйтйф. (20)

п

Заметим, что выражение в (20), стоящее под знаком интеграла, есть произведение функций, каждая из которых зависит лишь от одной переменной, что весьма существенно. Применять формулы типа (17) становится неудобно, так как граница 7 будет проходить внутри части

секторов (см. п. 5.3). Произведём универсализацию вычисления интегралов в (20). Для их вычисления применим формулу Грина-Стокса:

(р (11 )ds =11 (rot1,1 )dQ,

J Y J J П

где 1 = {P,Q, О}, P = P(х,у), Q = Q(x, у), ~1~ единичный вектор касательной к

кривой y, ограничивающей область Q, к - единичный вектор, перпендикулярный плоскости области Q. Линейная форма имеет вид:

С1,1 )ds = Pdx + Qdy = Pr dr + rQ^dp,

где Pr = P cosip + Q sin p, Qv = — P sinp + Q cos p. Выражение для ротора в полярной системе координат:

-V I l д . „ . l dPr м

rota = ~ — (rQv)-------------к

ror r op

Поэтому

cij = II Ci(p )Dj (r)rdrdp = JJ ^ “ dir (rQp) — ^ rdrdp = j) Pr dr + rQvdp.

Отсюда получаем, что

' l д , „ , l dPr

r or r o

Этому уравнению удовлетворяют

/і д l OP \

(-& (Q) — - = Ci(p)D (r). (21)

lr

Pr = О, Qv = - Ci(p) Dj (t)tdt.

r0

Иными словами, в качестве функции тQф(ф, т) возьмём первообразную от функции тБ/(т) (по т), умноженную на Сг(ф). Заметим, что эта первообразная есть сплайн, состоящий из п+2

выполнялось Qф(ф, 0) = 0. Отсюда получаем

c

ij

j) Ci(p) ^ J Dj(t)tdt^ dp. (22)

Обратим внимание на то, что непериодический фундаментальный сплайн Dj (т) = 0 при т < т^ если точка с координатами (|Pi,rj) не принадлежит некоторой области О$ D О. Поэтому Qф(ф,т) = 0 при т < т/. Итак, показано, что все коэффициенты еч равны нулю для таких пар (і, і), при которых точки с координатами ( фі, т/) £ О^ где 5 = 5(М, т, К).

5.7. Частный случай «простой» области

Область назовём «простой», если внутри неё найдётся такая точка, что любой луч, выпущенный из этой точки, пересечёт границу области только в одной точке. Поместим начало координат в эту точку и введём полярную систему координат. Тогда граница 7 области О задается функцией т = т(ф), ф Є [0, 2п] .

Зафиксируем некоторое ф. Заметим, что Dj(т) = 0 при т < тj—М2К.2, где М2-колпчество точек осреднения, используемых при построении Dj(т). Пусть < т/-М2Ъ,2.Тотда Dj(т) =

0 при т < Пусть т(ф) Є [£і2, £і2+і) (заметим, что І2 = І2(ф) зависит от угла ф и границы О

Гг(ф) І2_1 ГІв+1 гЛф')

/ tDj(і)М =^ / tDj(і)М + tDj(і)М =

]° з=іі]Ь К

І2 — 1 г Н2 5 Гг(ф)-&2 5

У2 / (и+^з)У2Ч3ия + (и+&2)У2з2и<1

3=11 7 0 д=0 д=0

12 1 сН2 гг(ф)-Іі2 _ 5

-'{3(и<1+1 + вН2ид)^ , . 7 ^ д

з=11'10 д=0 ■’0 д=0

" { Н2 , 1 \Ф)—‘ІІ2

Ьдз(ид+1 + зН2ид)йи + 'Е^Ьд2(ид+1 + І2Н2ид)йи

Е £ кзщ+2( -^ + _^\ + £ ці2((т(ф - І2Н2)++12Н2(т(ф) - І2Н2)+1 ^ г~1^ д=0 д 2 \9 + 2 9 + Ч 7=0 9 V 9 + 2

Поэтому

. 0 ч9 + 2 д + 1 ^ д \ д + 2 д + 1 у

з=11 д=0 4 1 1 ' д=0 4 '

°ч = І £ £ д^ + д+т) Сі(ф)^ф+

(т(ф) - І2 (ф)Н2)д+2 < „ (т(ф) - І2(ф)Н2)д+1\

^ 1=0

где

Ь1-1 5

Сг(ф) = £ Т.аре"'

п=0 р=0

Здесь арП и Ъ^ - р-й и д-й коэффициенты п-го и 8-го полиномов в г-м фундаментальном периодическом сплайне на окружности [0, 2п] и в И-м фундаментальном непериодическом сплайне на отрезке [0, К^. Шаг Н1 = 2п/Ь^ Фундаментальный периодический сплайн Сг(ф) строится по набору данных {ук = 5гкI к = 0,1,..., К\}. Шаг Н2 = 1/Ь2, фундаментальный непериодический сплайн Б/(т) строится по набору данных {ук = 5/кI к = 0,1,..., К2; у0,..., у^^1}, где у0,..., у0Р принимают значения либо 0, либо 1.

5.8. Оценка точности квадратурной формулы для двумерных односвязных областей

Обозначим через К = тах(Ь,1, Ь,2)- Пусть выполнены условия устойчивости матрицы и, например, т1/М1 < (,*, т2/М2 < и пусть / Є С(п+1)(Ог), где О$ Э О, т.е. мы предполагаем, что функция / определена и (п+1) раз непрерывно дифференцируема в несколько большей области О$ Э О. Поместим область О$ в круг К радиуса К. Введём полярную систему координат, взяв за начало координат центр круга К. Продолжим функцию / в К \ О$ тождественным нулем. Обозначим через Б(ф, т) т — ф-сплайн, приближающий таким обра-

/К

Теорема 5. Пусть Б(ф,т) - это т — ф-сплайн, приближающий функцию / , пусть

(М + т)К < р(^ё ,ч)- Здесь р(^ё, 7) ~ расстояние между границами областей О$ и О соответственно. Тогда справедлива оценка:

Кі — 1 К2

/(х,у)йО - Е Ё0"3Уі і=0 3=0

< ан(п+1\

(24)

п

Здесь уг/ = /(фг,т/) - значения функции / в узлах сетки, весовые коэффициенты сг/ определены формулами (22), (23), суммирование поизводится лишь по тем индексам г и И, для которых (фг ,т/) € П$.

Доказательство. Заметим, что

К1-1 К2

= 1

г=0 /=0

т.к. Б(ф,т) = 1, если / = 1. Из (13) следует, что \Б(ф,т) — /(ф,т)\ < С00Н(п+1'). Поэтому

К1-1 К2

/(х, у)йП — Е £ сг/уг/

г=0 /=0

<

/йП — БйП

пп

+

К1-1 К2

БйП —^ ^г/ уг/

г=0 /=0

<

< С00Н(п+1) шез(П) +

К1 - 1 К2

Цп (б — £ Т^уг/ Сг(фЩ (т) I йП

г=0 /=0

Последнее слагаемое равно нулю, так как Сг(ф)Б/(т) = 0 в области П для тех пар индексов ^ таторых (фг, т/) € П$. □

Замечание 1. Если заданную функцию / приближать ф — т-сплайном, то оценка (24) также будет справедлива, так как на круге К ф — т-сплайн отличается от т — ф-сплайна на величину 0(Н(п+1)).

Замечание 2. Особенно удобными получаются квадратурные формулы при использовании полулокальных сглаживающих сплайнов класса С0 при ш = 1.

Литература

1. Бабенко, К.И. Основы численного анализа / К.И. Бабенко. - М.; Ижевск: НИЦ Регулярная и хаотическая динамика, 2002.

2. Соболев, С.Л. Введение в теорию кубатурных формул / С.Л. Соболев. - М.: Наука, 1974.

3. Мысовских, И.П. Интерполяционные кубатурные формулы / М.П. Мысовских. - М.: Наука, 1981.

4. Крылов, А.Н. Лекции о приближенных вычислениях / А.Н. Крылов. - М.; Л.: ГИТТЛ, 1950.

5. Стечкин, С.Б. Сплайны в вычислительной математике / С.Б. Стечкин, Ю.Н. Субботин.

- М.: Наука, 1976.

6. Завьялов, Ю.С. Методы сплайн-функций / Ю.С. Завьялов, Б.И. Квасов, В.Л. Мирошниченко. - М.: Наука, 1980.

7. Колмогоров, А.Н. О представлении непрерывных функций нескольких переменных в виде суперпозиции функций одного переменного и сложения / А.Н. Колмогоров // Избранные труды. Математика и механика. - М.: Наука, 1985.

8. Соболев, С.Л. Кубатурные формулы / С.Л. Соболев, В.Л. Васкевич. - Новосибирск: Изд-во ИМ СО РАН, 1996.

9. Рамазанов, М.Д. Теория решетчатых кубатурных формул с ограниченным пограничным слоем / М.Д. Рамазанов. - Уфа: ИМВЦ УНЦ РАН, 2009.

п

п

10. Силаев, Д.А. О кубатурных формулах высокого порядка аппроксимации для широкого класса областей / Д.А. Силаев, Д.О. Коротаев // Сб. тр. XVI междупар. конф. «Математика. Компьютер. Образованно/ под ред. Г.Ю. Рпзнпченко. - Ижевск, 2009. - Т. 2.

- С. 20-38.

11. Силаев, Д.А. О кубатурных формулах высокого порядка аппроксимации для произвольных областей / Д.А. Силаев // Сб. тр. междун. конф. «Современная математика и математическое образование, проблемы истории и философии математики:»/ под ред. А.А. Артемова. - Тамбов, 2008. - С. 65-70.

Б

исследования и моделирование. - 2010. - Т. 2, № 4. - С. 349-358.

13. Силаев, Д.А. Приближение Б-сплайнамп гладких функций / Д.А. Силаев, Г.И. Якушина // Тр. семинара имени И.Г. Петровского. - 1984. - Вып. 10. - С. 197-206.

14. Силаев, Д.А. Б-сплайн на круге / Д.А. Силаев, Д.О. Коротаев // Тез. междупар. конф. «Математика. Компьютер. Образование». - Пущино, 2003. - С. 157.

Дмитрий Алексеевич Силаев, кандидат физико-математических наук, кафедра «Общие

»

ный университет им. М.В. Ломоносова, [email protected].

Bulletin of the South Ural State University. Series «Mathematical Modelling, Programming & Computer Software:»,

2013, vol. 6, no. 4, pp. 87-100.

MSC 65D32

Quadrature Formulas with High Order Approximation

D.A. Silaev, Lomonosov Moscow State University, Moscow, Russian Federation, [email protected]

In the article the method of creation the quadrature formulas with high order approximation for a wide class of the areas is given. This method is based on approach of smooth function on the plane by the semilocal smoothing spline or S-spline. Semilocal smoothing splines are initiated by D.A. Silaev. Earlier the splines of the third and fifth degree are considered and applied. This work is devoted to use of S-splines of higher degrees. Steady ^-splines of a class of C0 (only continuous), consisting of polynoms of high degree of n (n = 9,10) makes it possible to receive quadrature formulas of the 10th and 11th orders of approximation. It is supposed that integrand function belongs to Cp class (to p = 10,11) in a bigger area, than initial area on which integration is conducted. It is also supposed that the border of area is set parametrically that helps to consider area border with a fine precision. Similar approach is possible for the construction of cubature formulas.

Keywords: an approximation; a spline; integrals; quadrature formulas; numerical

methods.

References

1. Babenko K.I. Osnovy chislennogo analiza [Fundamentals of Numerical Analysis]. Moscow, Izhevsk, NITs Regulyarnaya i khaoticheskaya dinamika, 2002.

2. Sobolev S.L. Vvedenie v teoriyu kuhaturnykh formul [Introduction to the Theory of Cubature Formulas]. Moscow, Nauka, 1974.

3. Mysovskikh I.P. Interpolyatsionnye kubaturnye formuly [Java Applet Formula]. Moscow, Nauka, 1981.

4. Krylov A.N. Lektsii о priblizhennykh vychisleniyakh [Lectures on Approximate Calculations]. Moscow, Leningrad, GITTL, 1950.

5. Stechkin S.B., Subbotin Yu.N. Splayny v vychislitel’noy matematike [Splines in Computational Mathematics]. Moscow, Nauka, 1976.

6. Zav’yalov Yu.S., Kvasov B.I., Miroshnichenko V.L. Metody splayn-funktsiy [Methods of Spline Functions]. Moscow, Nauka, 1980.

7. Kolmogorov A.N. О predstavlenii nepreryvnykh funktsiy neskol’kikh peremennykh v vide superpozitsii funktsiy odnogo peremennogo i slozheniya [On the Representation of Continuous Functions of Several Variables by Superposition of Functions of One Variable and Addition]. Moscow, Nauka, 1985.

8. Sobolev S.L., V.L. Vaskevich V.L. Kubaturnye formuly [Cubature Formula]. Novosibirsk, Izd-vo IM SO RAN, 1996.

9. Ramazanov M.D. Teoriya reshetchatykh kubaturnykh formul s ogranichennym pogranichnym sloem [The Theory of Lattice Rules with a Limited Boundary Layer]. Ufa, IMVTs UNTs RAN, 2009.

10. Silaev D.A., Korotaev D.O. Cubature Formulas of High-Order Methods for a Wide Range of Areas [O kubaturnykh formulakh vysokogo poryadka approksimatsii dlya shirokogo klassa oblastey]. Sbornik trudov XVI mezhdunarodnoy konferentsii «Matematika. Komp’yuter. Obrazovanie» [Proceedings Works of the XVI International Conference flqq Mathematics. Computer. Education»], Izhevsk, 2009, vol. 2, pp. 20-38.

11. Silaev D.A. Cubature Formulas of High-Order Methods for Arbitrary Domains [O kubaturnykh formulakh vysokogo poryadka approksimatsii dlya proizvoPnykh oblastey].

«

»

«

»

12. Silaev D.A. Semilocal Smoothing S— splines [Polulokal’nye sglazhivayushchie S—splayny]. Komp’yuternye issledovaniya i modelirovanie [Computer Research and Modelling], 2010, vol. 2, no. 4, pp. 349-358.

13. Silaev D.A., Yakushina G.I. S-Spline Approximation of Smooth Functions [Priblizhenie S-splaynami gladkikh funktsiy]. Trudy seminara imeni I.G. Petrovskogo [Proceedings of the Seminar Named I.G. Petrovsky], 1984, issue 10, pp. 197-206.

14. Silaev D.A., Korotaev D.O. S-Spline Lap [S-splayn na kruge]. Tezisy mezhdunarodnoy

«»

«»

Поступила в редакцию 6 июня 2013 г.