Вероятностные методы. Метод Монте-Карло

Раздел № 9

В отличие от традиционных численных методов решения задач, заключающихся в разработке алгоритма построения решения, для исследования некоторых классов задач оказывается более целесообразным моделирование их сущности с использованием законов больших чисел теории вероятностей. Оценки f₁, f₂, ... , f_n искомой величины f в этом случае строят на основании статистической обработки материала, полученного в результате многократных случайных испытаний. Основным требованием при этом является сходимость по вероятности рассматриваемой случайной величины к искомому решению задачи:

(9.1)

Здесь Р – вероятность, ε – сколь угодно малое положительное число.

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

9.1. Вычисление определенного интеграла

Пусть задана непрерывная случайная величина ξ (кси), с плотностью вероятности р(х), значения которой распределены на интервале (а, b)

Плотность вероятности р(х) обладает следующими свойствами:

1) p(x) ³ 0;

2) . (9.2)

Тогда математическое ожидание случайной величины ξ равно интегралу

M(ξ) = . (χ – хи)

Для функции f(ξ), аргументом которой является случайная величина ξ, т. е. для случайной функции, математическое ожидание равно

M(f((ξ)) =. (9.3)

Отсюда следует, что любой интеграл вида

где функция р(х) обладает свойствами (9.2), можно принять за математическое ожидание некоторой случайной функции f(ξ).

Но математическое ожидание случайной величины f(ξ) можно приближенно определить с помощью статистических оценок, т. е. на основе выборки {ξ₁, ξ₂,..., ξ_N} объема N как среднее арифметическое значений {f(ξ₁), f(ξ₂),..., f(ξ_N)} поэтому интеграл (9.3) можно вычислить приближенно по формуле

» (9.4)

Теоретическим основанием для такого перехода является центральная предельная теорема теории вероятностей, которая в упрощенной формулировке утверждает следующее.

Среднее арифметическое N испытаний случайной величины ξ

ξ_N =

является случайной величиной с тем же математическим ожиданием

M(ξ_N) = М(ξ)

и дисперсией, равной

D(ξ_N) =

и при N ® ¥ закон распределения случайной величины ξ_N стремится к нормальному закону.

Очевидно, что чем больше N, тем меньше дисперсия D(ξ_N) = . Величину погрешности формулы (9.4)можно оценить по вероятности. Например, при больших N можно утверждать, что с вероятностью 0,997 ошибка не превосходит величину 3σ = 3×(правило «трех сигм» для нормально распределенной случайной величины).

Можно считать, что погрешность формулы (9.4) есть величина порядка O, но для повышения точности в данном случае нельзя применять правило Рунге-Ромберга.

Приведем другой способ статистической оценки для одномерного интеграла

.

Для этого вспомним его геометрический смысл. Предположим, что функция f(x) положительна на отрезке [a, b]. Тогда интеграл равен площади криволинейной трапеции, ограниченной графиком функции f(x), осью абсцисс и прямыми х = а и х = b (рис. 9.1).

Рассмотрим две случайные величины ξ – равномерно распределенную на отрезке [а, b] и η (эта) – равномерно распределенную на отрезке [0, f_max] где f_max = .

Рис. 9.1

Вероятность попадания случайной точки (ξ, η) в криволинейную трапецию равна отношению площади трапеции к площади прямоугольника

{(х, у), a £ х £ b, 0 £ у £ f_max }:

(9.5)

Проведем серию из N испытаний и получим N случайных точек (х, у), принадлежащих прямоугольной области {a £ х £ b, 0 £ у £ f_max }. Подсчитаем число N_f точек, удовлетворяющих условию у £ f(x). Тогда вероятность попадания случайной точки (ξ, η) в криволинейную трапецию приближенно равна относительной частоте попадания в криволинейную трапецию, т. е. р » и интеграл приближенно вычисляется по формуле

» (9.6)

Другой, более простой способ вычисления интеграла заключается в следующем.

Проведем серию из N испытаний случайной величины, равномерно распределенной на отрезке [а, b], и получим N случайных чисел x_i, принадлежащих отрезку [а, b]. Вычислим интеграл по формуле

» (9.7)

Отметим, что в (9.7) подынтегральная функция может принимать положительные и отрицательные значения, тогда как формула (9.6) применима только для неотрицательной функции f(x).

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

(9.8)

и применить формулу (9.7).