Моделирование случайных воздействий на системы

При моделировании системы 5 методом имитационного моделирования, в частности методом статистического моделирования на ЭВМ, существенное внимание уделяется учету случайных факторов и воздействий на систему. Для их формализации используются случайные события, дискретные и непрерывные величины, векторы, процессы. Формирование на ЭВМ реализаций случайных объектов любой природы из перечисленных сводится к генерации и преобразованию последовательностей случайных чисел. Вопросы генерации базовых последовательностей псевдослучайных чисел {х,}, имеющих равномерное распределение в интервале (0, 1), были рассмотрены в § 4.2, поэтому остановимся на вопросах преобразования последовательностей случайных чисел {х,} в последовательность {у,} для имитации воздействий на моделируемую систему 5.

Эти задачи очень важны в практике имитационного моделирования систем на ЭВМ, так как существенное количество операций, а значит, и временных ресурсов ЭВМ расходуется на действия со случайными числами. Таким образом, наличие эффективных методов, алгоритмов и программ формирования, необходимых для моделирования конкретных систем последовательностей случайных чисел {>>,}, во многом определяет возможности практического использования машинной имитации для исследования и проектирования систем [37, 46].

Моделирование случайных событий.

Простейшими случайными объектами при статистическом моделировании систем являются случайные события. Рассмотрим особенности их моделирования [4].

Пусть имеются случайные числа х_ь т. е. возможные значения случайной величины ?, равномерно распределенной в интервале (О,.

1). Необходимо реализовать случайное событие А, наступающее с заданной вероятностью р. Определим А как событие, состоящее в том, что выбранное значение х, случайной величины? удовлетворяет неравенству.

Тогда вероятность события А будет Р (А) = 1 4х=р. Противоположное событие А состоит в том, что х,>/?. Тогда Р (А)= 1 — р.

Процедура моделирования в этом случае состоит в выборе значений х_{ и сравнении их с р. При этом, если условие (4.16) выполняется, исходом испытания является событие А.

Таким же образом можно рассмотреть группу событий. Пусть А₁₉ А₂, А, — полная группа событий, наступающих с вероятностями /?!, /?₂, …, р_л соответственно. Определим А_т как событие, состоящее в том, что выбранное значение х, случайной величины {удовлетворяет неравенству.

где.

Тогда.

Процедура моделирования испытаний в этом случае состоит в последовательном сравнении случайных чисел x_i со значениями 4. Исходом испытания оказывается событие А_т, если выполняется условие (4.17). Эту процедуру называют определением исхода испытания по жребию в соответствии с вероятностями р_1У /*2″ •**" А*.

Эти процедуры моделирования были рассмотрены в предположении, что для испытаний применяются случайные числа х" имеющие равномерное распределение в интервале (0, 1). При моделировании на ЭВМ используются псевдослучайные числа с квазиравномерным распределением, что приводит к некоторой ошибке. Оценим ее Пример 4.7. Пусть имеются л-разрядные случайные числа с возможными значениями **" //(2*—1), /=0,2*-1. Подставив в (4.16) вместо х, число х], определим А* как событие, состоящее в том, что х* 4р- «

Вероятность наступления события А* может быть определена как Р (А*)*=т/2, где т — количество случайных чисел, меньших или равных р. Отсюда следует, что использование чиста х* вместо х_{ приводит к ошибке в определении вероятности события Ар^т^—р.

Очевидно, что максимальное значение ошибки не превосходит величины 1/(2л-1). Таким образом, для уменьшения влияния ошибок можно воспользоваться увеличением разрядности случайных чисел.

При моделировании систем часто необходимо осуществить такие испытания, при которых искомый результат является сложным событием, зависящим от двух (и более) простых событий. Пусть, например, независимые события А и В имеют вероятности наступления р_л и р_в. Возможными исходами совместных испытаний в этом случае будут события АВ, ИВ, АВ, ИВ с вероятностями р^_ъ, О —Рл)Рв. Рл ( 1 -Рв), (1 -Рл)( 1 -Рв)•

Для моделирования совместных испытаний можно использовать два варианта процедуры: 1) последовательную проверку условия.

• (4.16); 2) определение одного из исходов АВ, ИВ, АВ, ИВ по жребию с соответствующими вероятностями, т. е. аналогия (4.17). Первый вариант требует двух чисел х₁ и сравнений для проверки условия
• (4.16). При втором варианте можно обойтись одним числом х₍, но сравнений может потребоваться больше. С точки зрения удобства построения моделирующего алгоритма и экономии количества операций и памяти ЭВМ более предпочтителен первый вариант.

Рассмотрим теперь случай, когда события А и В являются зависимыми и наступают с вероятностями р_А и р_в. Обозначим через Р (В/А) условную вероятность наступления события В при условии, что событие А произошло. При этом считаем, что условная вероятность Р (В/А) задана.

Рассмотрим один из вариантов построения модели. Из последовательности случайных чисел {*,} извлекается очередное число х_т и проверяется справедливость неравенства х_т<�р_А. Если это неравенство справедливо, то наступило событие А. Для испытания, связанного с событием В, используется вероятность Р (В/А). Из совокупности чисел {*,} берется очередное число х_т+_х и проверяется условие х_т+_{^Р (В/А). В зависимости от того, выполняется или нет это неравенство, исходом испытания являются А В или АВ.

Если неравенство х_т<�р_А не выполняется, то наступило событие А. Поэтому для испытания, связанного с событием В, необходимо определить вероятность.

Выберем из совокупности {*,} число х_т+: и проверим справедливость неравенства х_т+1 ^ Р (В[Л). В зависимости от того, выполняется оно или нет, получим исходы испытания АВ или А В.

Логическая схема алгоритма для реализации этого варианта модели показана на рис. 4.13. Здесь ВИД[…] — процедура ввода исходных данных; ГЕН […] — генератор равномерно распределенных случайных чисел; ХМ=х_т ХМЫх_т^; РА=р_л; РВшр_ш РВА = Р (В/А)', РВЫА = Р (В/А)] КА, КЯА_У КАВ_У КАЫВ, КАТАВ, ЮТАМ В — число событий А, Л_у АВ, ХВ, ~АВ_у А В соответственно; ВРМ […] — процедура выдачи результатов моделирования.

Рассмотрим особенности моделирования на ЭВМ марковских цепей, служащих, например, для формализации.

Рис. 4.13. Схема моделирующего алгоритма при зависимых событиях.

процессов в Р-схемах (см. § 2.4).

Простая однородная марковская цепь определяется матрицей переходов.

где рц — вероятность перехода из состояния z_i в состояние Матрица переходов Р полностью описывает марковский процесс. Такая матрица является стохастической, т. е. сумма элементов каждой строки равна единице:

Обозначим через р,(п)_у 1=1, к, вероятности того, что система будет находиться в состоянии г, после п переходов. По определению.

Используя событийный подход, можно подойти к моделированию марковской цепи следующим образом. Пусть возможными исходами испытаний являются события А_и Л₂, А*. Вероятность Р/у — это условная вероятность наступления события Ау в данном испытании при условии, что исходом предыдущего испытания было событие А{. Моделирование такой цепи Маркова состоит в последовательном выборе событий Aj по жребию с вероятностями рф

Сначала выбирается начальное состояние г*, задаваемое начальными вероятностями рЛ0), р₂(0), •••" Р*(0). Для этого из последовательности чисел {*,} выбирается число х_т и сравнивается с /, из (4.17), где в качестве /?, используются значения р₁ (0), р₂(0)> …, р_к(0). Таким образом, выбирается номер т₀, для которого оказывается справедливым неравенство (4.17). Тогда начальным событием данной реализации цепи будет событие А«₀. Затем выбирается следующее случайное число х_{т+ ь} которое сравнивается с /_г, где в качестве р₁ используются р_то_;. Определяется номер т_1? и следующим событием данной реализации цепи будет событие А_т и т. д. Очевидно, что каждый номер т, определяет не только очередное событие А_ш формируемой реализации, но и распределение вероятностей р_т1и р"12з…» Рт* для выбора очередного номера т_/+1, причем для эргодических марковских цепей влияние начальных вероятностей быстро уменьшается с ростом номера испытаний. Эргодическим называется всякий марковский процесс, для которого предельное распределение вероятностей р, (л), /= 1, к, не зависит от начальных условий р_/(0). Поэтому при моделировании можно принимать, что.

Аналогично можно построить и более сложные алгоритмы, например, для моделирования неоднородных марковских цепей [29].

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

Для формирования возможных значений случайных величин с заданным законом распределения исходным материалом служат базовые последовательности случайных чисел {*,}, имеющие равномерное распределение в интервале (0, 1). Другими словами, случайные числа как возможные значения случайной величины, имеющей равномерное распределение в интервале (0, 1), могут быть преобразованы в возможные значения у_] случайной величины гI закон распределения которой задан [4].

Моделирование дискретных случайных величии. Рассмотрим особенности преобразования для случая получения дискретных случайных величин. Дискретная случайная величина у принимает значения у₁ с вероятностями р₁₉ р_г>…, ру, составляющими дифференциальное распределение вероятностей.

При этом интегральная функция распределения.

Для получения дискретных случайных величин можно использовать метод обратной функции [10, 29, 53]. Если { — равномерно распределенная на интервале (0, 1) случайная величина, то искомая случайная величина г/ получается с помощью преобразования.

где ¹ — функция, обратная Е_п.

Алгоритм вычисления по (4.19) и (4.20) сводится к выполнению следующих действий:

При счете по (4.21) среднее число циклов сравнения.

Пример 4.8. Необходимо методом обратной функции на основании базовой последовательности случайных чисел {х;}, равномерно распределенных в интервале (0, 1), получить последовательность чисел (уу), имеющих биномиальное распределение, задающее вероятность у удачных исходов в N реализациях некоторого эксперимента:

где р=0,5 и N=6; С*" ЛГ!/у!(Л-у)!

Математическое ожидание и дисперсия биномиального распределения соответственно будут Му] = пр, Ву] = пр (1 -р).

Используя для ру обозначения, принятые" (4.21), вычислим:

Например, получив из равномерного распределения число *,=0,85 393 и проведя сравнения по алгоритму (4.21), найдем, что х/*0,85 393 0,89 063, т. е. уу•*4.

При этом среднее число циклов сравнения ?*" 1 0,1 562+2 0,9 375 ч- + 3 0,23 438 + 4 0,31 250 + 5 0,23 438 + 6 (0,9 375 4−0,1 562)" 3,98.

Пример 4.9. Необходимо проверить стохастичность последовательности из N случайных чисел {уу}, полученных при имитации биномиального распределения при заданных параметрах пир. Простейшим способом проверки является оценка выполнения условий.

Проверим на соответствие биномиальному распределению с параметрами п «5 и р=0,1 такой последовательности случайных чисел: (уу} = 0, 0, 1,0, 1, 0, 2, 0, 1,0; ЛГ- 10.

Вычислим л/? «0,5; пр ( 1 4-пр-р) —0,7;

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

Можно привести и другие примеры алгоритмов и программ получения дискретных случайных величин с заданным законом распределения, которые находят применение в практике моделирования систем на ЭВМ.

Моделирование непрерывных случайных величии. Рассмотрим особенности генерации на ЭВМ непрерывных случайных величин. Непрерывная случайная величина г} задана интегральной функцией распределения.

где/_п{у) — плотность вероятностей.

Для получения непрерывных случайных величин с заданным законом распределения, как и для дискретных величин, можно воспользоваться методом обратной функции. Взаимно однозначная монотонная функция ^ = полученная решением относительно г уравнения Р"(у) = ?, преобразует равномерно распределенную на интервале (0, 1) величину { в г с требуемой плотностью/_п(у).

Действительно, если случайная величина г имеет плотность распределения /_ч(у), то распределение случайной величины.

является равномерным в интервале (0, 1) [4, 29]. На основании этого можно сделать следующий вывод. Чтобы получить число, принадлежащее последовательности случайных чисел {у^, имеющих функцию плотности /,(у), необходимо разрешить относительно у₎ уравнение.

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

Пример 4.10. Необходимо получить случайные числа с показательным законом распределения.

В силу соотношения (4.22) получим.

где х₁ — случайное число, имеющее равномерное распределение в интервале (0, 1). Тогда.

Учитывая, что случайная величина? имеет также равномерное распределение в интервале (0, 1), можно записать уу* — ^ In x_t.

Пример 4.11. Необходимо получить случайные числа у_у с законом распределения.

Воспользовавшись (4.22), получим X (у,—Ду?У4)=х,. Отсюда.

Можно привести и другие примеры использования соотношения (4.22). Но этот способ получения случайных чисел с заданным законом распределения имеет ограниченную сферу применения в практике моделирования систем на ЭВМ, что объясняется следующими обстоятельствами: 1) для многих законов распределения, встречающихся в практических задачах моделирования, интеграл (4.22) не берется, т. е. приходится прибегать к численным методам решения, что увеличивает затраты машинного времени на получение каждого случайного числа; 2) даже для случаев, когда интеграл (4.22) берется в конечном виде, получаются формулы, содержащие действия логарифмирования, извлечения корня и т. д., которые выполняются с помощью стандартных подпрограмм ЭВМ, содержащих много исходных операций (сложения, умножения и т. п.), что также резко увеличивает затраты машинного времени на получение каждого случайного числа.

Поэтому в практике моделирования систем часто пользуются приближенными способами преобразования случайных чисел, которые можно классифицировать следующим образом: а) универсальные способы, с помощью которых можно получать случайные числа с законом распределения любого вида; б) не универсальные способы, пригодные для получения случайных чисел с конкретным законом распределения [4, 46].

Рассмотрим приближенный универсальный способ получения случайных чисел, основанный на кусочной аппроксимации функции плотности. Пусть требуется получить последовательность случайных чисел да с функцией плотности /_ч(у), возможные значения которой лежат в интервале (а, Ь). Представим /_ч(у) в виде кусочно-постоянной функции, т. е. разобьем интервал (а, Ь) на т интервалов, как это показано на рис. 4.14, и будем считать /^(у) на каждом интервале постоянной. Тогда случайную величину г] можно представить в виде Т} = а_к+г!_к*, где а_к — абсцисса левой границы ?-го интервала; т_к* — случайная величина, возможные значения которой располагаются равномерно внутри к-го интервала, т. е. на каждом участке а*ч-а*₊₁ величина г]_к* считается распределенной равномерно. Чтобы аппроксимировать/,(у) наиболее удобным для практических целей способом, целесообразно разбить (а, Ь) на интервалы так, чтобы вероятность попадания случайной величины т/ в любой интервал (а_ь а_к+_х) была постоянной, т. е. не зависела от номера интервала к. Таким образом, для вычисления а_к воспользуемся следующим соотношением:

Алгоритм машинной реализации этого способа получения случайных чисел сводится к последовательному выполнению следующих действий: 1) генерируется случайное равномерно распределенное число х₍ из интервала (О, '1); 2) с помощью этого числа.

Рис. 4.14. Кусочная аппроксимация функции плотности.

случайным образом выбирается интервал (д*, 1); 3) генерируется число *1+1 и масштабируется с целью приведения его к интервалу (а_к, а_к+,), т. е. домножается на коэффициент (д*+1 — Д*)*/*. ь 4) вычисляется случайное число у^=а_к+(а_к+х—а_к)дг_/+1 с требуемым законом распределения.

Рассмотрим более подробно процесс выборки интервала (а_к, й_к+) с помощью случайного числа х,. Целесообразно для этой цели построить таблицу (сформировать массив), в которую предварительно поместить номера интервалов к и значения коэффициента масштабирования, определенные из соотношения (4.23) для приведения числа к интервалу (а, Ь). Получив из генератора случайное число х₁₉ с помощью этой таблицы сразу определяем абсциссу левой границы а_к и коэффициент масштабирования (д*+| — а_к).

Достоинства этого приближенного способа преобразования случайных чисел: при реализации на ЭВМ требуется небольшое количество операций для получения каждого случайного числа, так как операция масштабирования (4.23) выполняется только один раз перед моделированием, и количество операций не зависит от точности аппроксимации, т. е. от количества интервалов т.

Рассмотрим способы преобразования последовательности равномерно распределенных случайных чисел {х,} в последовательность с заданным законом распределения {#} на основе предельных теорем теории вероятностей. Такие способы ориентированы на получение последовательностей чисел с конкретным законом распределения, т. е. не являются универсальными [29, 43]. Поясним сказанное примерами.

Пример 4.12. Пусть требуется получить последовательность случайных чисел {у,}, имеющих нормальное распределение с математическим ожиданием а и средним квадратическим отклонением <�г.

Будем формировать случайные числа уу в виде сумм последовательностей случайных чисел {*/}, имеющих равномерное распределение в интервале (0, 1). Так как исходной (базовой) последовательностью случайных чисел {х<} при суммировании является последовательность чисел, имеющих равномерное распределение в интервале (0, 1), то можно воспользоваться центральной предельной теоремой для одинаково распределенных случайных величин (см. § 4.1): если независимые одинаково распределенные случайные величины х,…, х" имеют каждая математическое ожидание а, и среднее квадратическое отклонение <�г_и то сумма? х_/ асимптотически нормальна с математическим ожиданием, а = ;Уа, и средним квадратическим отклонением а = а и N. м

Как показывают расчеты, сумма? х, имеет распределение, близкое к нормальном}', уже при сравнительно небольших N. Практически для получения последовательности нормально распределенных случайных чисел можно пользоваться значениями #="8 — 12, а в простейших случаях — меньшими значениями ЛР, например, ЛГ=4+5 [4].

Пример 4.13. Пусть необходимо получить случайные числа, имеющие закон распределения Пуассона:

Для этого воспользуемся предельной теоремой Пуассона (см. § 4.1): если Р — вероятность наступления события А при одном испытании, то вероятность наступления т событий в N независимых испытаниях при N-* оо, p~*Q, Np = X асимптотически равна Р (т).

Выберем достаточно большое N, такое, чтобы p*=X/N< 1, будем проводить серии по N независимых испытаний, в каждом из которых событие А происходит с вероятностью р, и будем подсчитывать число случаев yj фактического наступления события А в серии с номером j. Числа yj будут приближенно следовать закону Пуассона, причем тем точнее, чем больше N. Практически N должно выбираться таким образом, чтобы р<0,1 +0,2 [4].

Алгоритм генерации последовательности случайных чисел ур имеющих пуассоновское распределение, с использованием данного способа показан на рис. 4.15. Здесь LA «Я; N~N PNsp; XImx_t — случайные числа последовательности, равномерно распределенной в интервале (0, 1); YJ^yj, NO — вспомогательная переменная; ВИД (…) — процедура ввода исходных данных; ВЫЧ (…) — процедура вычисления; ГЕН (…) — процедура генерации случайных чисел; ВРМ (…) — процедура выдачи результатов моделирования.