Снова рассмотрим интеграл по K-мерной области, разбитой сеткой на ячейки (рис. 2). Его можно вычислить последовательным интегрированием:

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

Последовательное интегрирование по всем направлениям приводит к кубатурным формулам, которые являются прямым произведением одномерных квадратурных формул:

(11)

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

Вообще говоря, для разных направлений можно использовать квадратурные формулы разных порядков точности . Тогда главный член погрешности имеет вид:

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

Можно подобрать веса и положение линий сетки так, чтобы одномерная квадратурная формула была точна для многочлена максимальной степени, т.е. была бы формулой Гаусса, тогда, для случая K=2:

(12)

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

Размещено на http://www.allbest.ru/

Произвольная область. Метод последовательного интегрирования можно применять к области произвольной формы, например, с криволинейной границей. Рассмотрим этот случай при K=2. Для этого проведём через область хорды, параллельные оси , и на них введём узлы, расположенные на каждой хорде так, как нам требуется (рис. 4). Представим интеграл в виде:

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

При вычислении интеграла по имеется одна тонкость. Если область ограничена гладкой кривой, то при длина хорды стремится к нулю не линейно, а как ; значит, вблизи этой точки . То же будет при . Поэтому интегрировать непосредственно по формулам высокого порядка точности бессмысленно. Целесообразно выделить из основную особенность в виде веса , которому соответствуют ортогональные многочлены Чебышева второго рода. Тогда второе интегрирование выполняется по формулам Гаусса-Кристоффеля:

(13)

где , а и -нули и веса многочленов Чебышева второго рода.

Чтобы можно было применять эту формулу, надо ординаты хорд на рис. 4 заранее выбрать в соответствии с узлами (13). Если это не было сделано, то придётся ограничиться интегрированием по обобщённой формуле трапеций, причём её эффективный порядок точности в этом случае будет ниже второго.

3. Метод Симпсона для кратных интегралов

Размещено на http://www.allbest.ru/

Пусть сначала область интегрирования есть K-мерный пространственный параллелепипед (рис. 5), стороны которого параллельны осям координат. Каждый из промежутков разобьём пополам точками:

,

где .

Всего таким образом, получим точек сетки.

Имеем:

кубатурный интеграл решение программа

(14)

Находим K-мерный интеграл, вычисляя каждый внутренний интеграл по квадратурной формуле Симпсона на соответствующем отрезке. Проведём полностью все вычисления для случая K=2:

Применяя к каждому интегралу снова формулу Симпсона, получим:

или

(15)

Формулу (15) будем называть кубатурной формулой Симпсона. Следовательно,

(15)

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

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

Размещено на http://www.allbest.ru/

Опять рассмотрим случай K=2. Положим, что стороны прямоугольника мы разделили соответственно на и равных частей; в результате получилась

относительно крупная сеть прямоугольников (на рис. 6 вершины этих прямоугольников отмечены более крупными кружками). Каждый из этих прямоугольников в свою очередь разделим на четыре равные части. Вершины этой последней мелкой сети прямоугольников примем за узлы кубатурной формулы.

Пусть и . Тогда сеть узлов будет иметь следующие координаты:

Для сокращения введём обозначение

Применяя формулу (15) к каждому из прямоугольников крупной сети, будем иметь (рис.6):

Отсюда, делая приведение подобных членов, окончательно находим:

(16)

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

Размещено на http://www.allbest.ru/

Если область интегрирования - произвольная, то строим параллелепипед , стороны которого параллельны осям координат (рис. 83). Рассмотрим вспомогательную функцию

В таком случае, очевидно, имеем:

Последний интеграл приближённо может быть вычислен по общей кубатурной формуле (16).

4. Практическая реализация метода

4.1 Пакетная реализация