Диссертация (1150625), страница 11
Текст из файла (страница 11)
Во-первых, оно должно быть совместимо с определением (,)последовательностей (определение 4) в том смысле, что если для фиксированного каждое измножеств разбиения будет состоять из элементарных множеств, то это гарантирует попаданиеодинакового количества точек в каждое из X1 , . . . , X . Во-вторых, алгоритм разбиения долженбыть таким, чтобы соблюдалась „вложенность“ ( , 2 )-разбиения в ( , 2+1 )-разбиение, тогда будет справедливо утверждение теоремы 10. В-третьих, разбиение должно быть достаточнопростым с алгоритмической точки зрения, то есть таким, чтобы задача определения индексаподмножества для произвольной точки ∈ была бы не очень вычислительно затратной.В работе автора и Ермакова С.М. [18], которая содержит первые численные результаты работы алгоритма Qint, предложено разбиение, определяемое только первой координатой .
Такое разбиение удовлетворяет всем вышеперечисленным требованиям, но, как было выясненопозднее, не является оптимальным. В самом деле, X1 , . . . , X представляют собой гиперпараллелепипеды, у которых с ростом уменьшается только одна грань. Такой алгоритм разбиения,по-видимому, хорошо подходит для функций, которые сильно изменяются только по первой координате.
Более универсальным оказывается алгоритм разбиения, предложенный автором в [19],именно его мы будем использовать во всех последующих численных экспериментах.Предлагаемый метод построения ( , 2 )-разбиений гиперкуба мы будем называть правилом бинарного рассечения. Определим шаги этого правила по индукции. Для = 0 ситуация тривиальна: X1 = задаёт ( , 1)-разбиение.
Для = 1 рассечём гиперкуб плоскостью, задаваемой уравнением 1 = 0.5, и тогда ( , 2)-разбиение будет состоять из множествX1 = [0, 0.5] × [0, 1] × . . . × [0, 1], X2 = [0.5, 1] × [0, 1] × . . . × [0, 1]. Теперь для = 2, . . . , 59(a) Шаг 1(b) Шаг 2(c) Шаг 3(d) Шаг 4(e) Шаг 5(f) Шаг 6Рисунок 2.1: Первые шесть шагов правила бинарного рассечения, = 3мы можем построить каждое последующее разбиение из предыдущего, проведя дополнительноерассечение плоскостью = 0.5.Опишем переход от ( , − 1)-разбиения к ( , )-разбиению для > .1) Если имеет остаток 1 от деления на (то есть представимо в виде + 1, ∈ N0 ),мы смотрим на рассечение, проведённое шагов назад. Оно проведено по первой координате несколькими плоскостями: 1 = 1 , 1 = 2 , .
. . , 1 = . Пусть числа 1 , 2 , . . . , упорядочены по возрастанию. Новое рассечение зададим плоскостями 1 =1 , 1 =2−1, 12= 2 , . . . , 1 = , 1 =1−.21, 12=Иными словами, каждый из отрезков[0, 1 ], [1 , 2 ], . . . , [ , 1] получает дополнительное рассечение пополам.2) В противном случае применяем рассечение для шага − 1, но по следующей координате.Номер этой координаты равен остатку от деления на (за тем лишь исключением, когдаон равен 0, тогда рассечение проводится по последней координате, ).Правило бинарного рассечения для случая = 3 проиллюстрировано рисунком 2.1, где изображены первые 6 шагов.Как уже было отмечено, оценивание дисперсии для рандомизированного квази-Монте-Карлопроводится на основе некоторого числа повторений. Каждое такое повторение представляет со-60бой один и тот же квазислучайный набор, рандомизированный по-новому.
Это необходимо, поскольку для одного повтора сколько-нибудь надежная процедура оценки дисперсии неизвестна,и в этом смысле необходимость проводить повторы является несомненным недостатком квазиМонте-Карло. Отметим, что и наивный Монте-Карло, и Qint такого недостатка лишены; болеетого, Qint задумывался как способ борьбы именно с необходимостью „лишних“ повторений.Вместе с тем, Qint можно использовать в сочетании с повторами, аналогичными квази-МонтеКарло.
На данном этапе нам кажется, что сравнение с квази-Монте-Карло необходимо проводить именно так, повторяя процедуру Qint столько же раз, сколько повторяется квази-МонтеКарло. Это количество мы будем именовать числом внешних повторов. Таким образом, мы будем сравнивать Qint с наивным Монте-Карло без внешних повторов, а с рандомизированнымквази-Монте-Карло – с внешними повторами, число которых будет равно 16 (этого достаточнов соответствии с практическими рекомендациями, указанными Лемьё К. [31]).
Дополнительноотметим, что вопрос о количестве внешних повторов и их целесообразности в целом не разреˆˆ обладают устойчивостью и при меньшемшён окончательно: как показывает практика, оценки количестве внешних повторов, что является ещё одним потенциальным улучшением метода Qint.2.5Результаты численных экспериментовВсе численные результаты будут представлены в виде графиков, где по оси отложено общееколичество вычислений подынтегральной функции, а по оси – оценка стандартного отклонения, полученная тем или иным рассматриваемым методом (наивный Монте-Карло, рандомизированный квази-Монте-Карло или Qint). Сразу стоит отметить, что все сравниваемые схемыявляются несмещёнными, поэтому нас не будет интересовать среднее (его математическое ожидание равно искомому значению интеграла, ; забегая вперёд, во всех экспериментах для всехметодов значение попадает в доверительный интервал с уровнем доверия 99%), а только стандартное отклонение.Для каждого из методов в этих координатах будем называть его траекторией последовательно соединённые точки, отвечающие разному количеству вычислений функции.
Траекторииодного и того же метода будут обозначаться одним цветом. Обе оси логарифмические, поэтомуможно легко определять скорость убывания стандартного отклонения для той или иной траек(︀ )︀тории, и, соответственно, метода: её можно оценить как 1 , если траектория параллельнапрямой вида = − + . Построение такого рода профилей является стандартной процедурой,позволяющей визуально оценивать скорость сходимости различных алгоритмов (см., например,Оуэн А. [45]).Для метода Qint мы будем рассматривать семейство траекторий вида „Qint_q“. Для каждойиз таких траекторий фиксирован параметр (он связан с , числом разбиений), а с увеличениемколичества вычислений функции наращивается параметр (он связан с , числом точек в каждом разбиении).
Минимально возможный параметр равен нулю: это означает, что разбиениевырождается в случай X1 = . Таким образом, процедура оценивания дисперсии в точности61совпадает с процедурой Монте-Карло, поэтому мы ожидаем, что траектория Монте-Карло итраектория „Qint_0“ будут совпадать.Наконец, мы будем выделять траекторию „Best Qint“. Она, напротив, соответствует фиксированному и увеличивающемуся . Более того, значение = 2 фиксировано для всех последующих функций и вариаций экспериментов.
Это значит, что в каждое из подмножеств X1 , . . . , Xпопадает по 4 точки. Такой выбор обеспечивает достаточное усреднение для оценивания величинˆˆ } . Мы могли бы рассмотреть и другие фиксированные , но нам достаточно траектории{=1„Best Qint“ для оценки скорости убывания дисперсии метода.Все результаты последующих численных экспериментов получены автором с использованием языка C++ [58] и открытой библиотеки HIntLib за авторством Шюрера Р. [42].
Обработкарезультатов и получение иллюстраций проводилось с использованием языка R [59] и ряда пакетов дополнений, в особенности randtoolbox [60] и ggplot2 [61]. Код программы, при помощикоторой проведены расчёты, распространяется автором диссертации под свободной лицензиейGPL (старше версии 2) и доступен в репозитории на сервисе GitHub [62].2.5.1Произведение кубических полиномовДля построения первого примера рассмотрим простейший кубический полином 3 + 0.75.В качестве подынтегральной функции будем использовать произведение таких полиномов покаждой из размерностей:1 () =∏︁(3 + 0.75).(2.35)=1Значение интеграла от функции 1 не зависит от размерности и равно единице:∫︁1 =1 () =∏︁=1(︃43+44)︃⃒1⃒⃒⃒⃒= 1.(2.36) =0Поведение стандартного отклонения квадратуры Qint для функции 1 в одномерном случаеописано нами в предыдущей главе.
Так, соотношение (2.31), полученное через оценивание модуля непрерывности (или, что равнозначно в плане порядка, интегрального модуля непрерывности(︁)︁для = 2), предполагает, что скорость убывания не может быть медленнее, чем 1 3/2.Результат численного эксперимента для функции 1 в одномерном случае приведён на иллюстрациях 2.2 и 2.3 (без внешних повторов и с внешними повторами, соответственно).Отметим несколько важных особенностей полученного результата.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
