книги / Моделирование систем

чивает возможности использования метода серединных квадратов.

Конгруэнтные процедуры генерации. Широкое применение при моделировании систем на ЭВМ получили конгруэнтные процедуры генерации псевдослучайных последовательностей, представляющие собой арифметические операции, в основе которых лежит фунда­ ментальное понятие конгруэнтности. Два целых числа а и /? конгру­ энтны (сравнимы) по модулю т, где т — целое число, тогда и толь­ ко тогда, когда существует такое целое число к, что a,—f}=km, т. е. если разность a — делится на т и если числа а и /?дают одинаковые остатки от деления на абсолютную величину числа т. Например, 1984= 4 (mod 10), 5008 s 8 (mod 103) и т. д.

Конгруэнтные процедуры являются чисто детерминированны­ ми, так как описываются в виде рекуррентного соотношения, когда функция (4.9) имеет вид

121

где Xh X, p, M — неотрицательные целые числа. Раскроем рекуррентное соотношение (4.10):

Хх= ХХ0+ р (mod М);

X2=XXL+ p = Х2Х 0+ (Я+1 )ju(mod М);

X3=XX2+p=X3X0+(l2+M )p=X3X0+(X3-l)p/(X-\)(modM);

Х,=i % + (Л1— 1М Л - l)(mod М).

(4.11)

Если заданы начальное значение Х09 множитель X н аддитивная константа р, то (4.11) однозначно определяет последовательность целых чисел {Х{}, составленную из остатков от деления на М членов

последовательности {Х*Х0+ц(Х*—1)/(Я—1)}. Таким образом, для любого 1 справедливо неравенство Xt<M . По целым числам последовательности {Х^\ можно построить последовательность

{*/}= {XJM} рациональных чисел из единичного интервала (0, 1). Конгруэнтная процедура получения последовательностей псев­

дослучайных квазиравномерно распределенных чисел может быть реализована мультипликативным либо смешанным методом.

Мультипликативный метод. Задает последовательность неотри­ цательных целых чисел {Xi}, не превосходящих М, по формуле

т.е. это частный случай соотношения (4.10) при ц = 0.

Всилу детерминированности метода получаются воспроизводи­ мые последовательности. Требуемый объем машинной памяти при

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

Для машинной реализации наиболее удобна версия M = pg, где р — число цифр в системе счисления, принятой в ЭВМ (р=2 для двоичной и р = 10 для десятичной машины); g — число битов в ма­ шинном слове. Тогда вычисление остатка от деления на М сводится к выделению g младших разрядов делимого, а преобразование целого числа Xt в рациональную дробь из интервала х{ е (0, 1) осуществляется подстановкой слева от Х {двоичной или десятичной запятой.

Алгоритм построения последовательности для двоичной маши­ ны М = 2в сводится к выполнению таких операций [31, 36, 46]:

122

1.Выбрать в качестве Х 0 произвольное нечетное число.

2.Вычислить коэффициент Л=8*+3, где t — любое целое поло­ жительное число.

3.Найти произведение ЛХ0, содержащее не более 2g значащих разрядов.

4.Взять g младпшх разрядов в качестве первого члена последо­ вательности Х 19 а остальные отбросить.

5.Определить дробь x 1= X J2s из интервала (0,1).

6.Присвоить X 0—X v

7.Вернуться к п. 3.

Пример 4.4. Необходимо получить числа последовательности для случая $=4, используя приведенный алгоритм мультипликативного метода. Для этого выполня­ ем следующие действия: 1. Выбираем ДТою—? (в десятичной системе счисления) или J o=0111 (в двоичной системе счисления). 2. Найдем *=1, тогда Д10=11 или 5; пусть ^0 = 5 , X«0101. 3. Рассчитываем произведение ХХ0, берем g младпшх разрядов, вычисляем Х х и присваиваем Х 0= Х х, т. е. выполняем п. 3 — 7 алгоритма:

а) AY0=(0101)(0111)—00100011, *,=0011, *.=3/16=0,1875; б) 2*!=(0101) (0011)=00001 111, * 2= 1111, х ,= 15/16=0,9375; в) АЗГ2=(0101)(Ш 1)=01001011, * 3=1011, х3=11/16=0,6875; г) А*3=(0101) (1011)=00110111, *4=0111, х4= 7/16=0,4375.

Смешанный метод. Позволяет вычислить последовательность неотрицательных целых чисел {Afj}, не превосходящих М 9 по формуле

Xi+1 =в ХХг+ д (mod М),

т. е. в отличие от мультипликативного метода /i# 0. С вычислитель­ ной точки зрения смешанный метод генерации сложнее мультип­ ликативного на одну операцию сложения, но при этом возможность выбора дополнительного параметра позволяет уменьшить возмож­ ную корреляцию получаемых чисел.

Качество конкретной версии такого генератора можно оценить только с помощью соответствующего машинного эксперимента.

В настоящее время почти все пакеты прикладных программ универсальных ЭВМ для вычисления последовательностей равно­ мерно распределенных случайных чисел основаны на конгруэнтной процедуре [17, 38, 40].

4.3. ПРОВЕРКА И УЛУЧШЕНИЕ КАЧЕСТВА ПОСЛЕДОВАТЕЛЬНОСТЕЙ ПСЕВДОСЛУЧАЙНЫХ ЧИСЕЛ

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

123

последовательность псевдослучайных чисел удовлетворяет предъяв­ ляемым к ней требованиям, так как в противном случае даже при наличии абсолютно правильного алгоритма моделирования процес­ са функционирования системы S по результатам моделирования нельзя достоверно судить о характеристиках системы [29, 37].

Проверка качества последовательностей. Результаты анализа си­ стемы S, полученные методом статистического моделирования на ЭВМ, существенно зависят от качества используемых псевдослучай­ ных квазиравномерных последовательностей чисел. Поэтому все применяемые генераторы случайных чисел должны перед моделиро­ ванием системы пройти тщательное предварительное тестирование, которое представляет собой комплекс проверок по различным ста­ тистическим критериям, включая в качестве основных проверки (тесты) на равномерность, стохастичность и независимость. Рас­ смотрим возможные методы проведения таких проверок, наиболее часто используемые в практике статистического моделирования систем.

Проверка равномерности последовательностей псевдослучайных квазиравномер­ но распределенных чисел {х,} может быть выполнена по гистограмме с использова­ нием косвенных признаков [4, 26]. Суть проверки по гистограмме сводится к следу­ ющему. Выдвигается гипотеза о равномерности распределения чисел в интервале (О, 1). Затем интервал (0, 1) разбивается на т равных частей, тогда при генерации

последовательности {х,} каждое из чисел х с вероятностью Pj= ljm ,j= 1, т, попадает

в один из подынтервалов. Всего в каждый j-й подынтервал попадает Nj чисел

____ т

последовательности {х,}, i= l, N, причем N= £ Nj. Относительная частота попада-

j- 1

ния случайных чисел последовательности {х,-} в каждый из подынтервалов будет равна Nj/N. Вид соответствующей гистограммы для примера показан на рис. 4.11, а, где пунктирная линия соответствует теоретическому значению pj, а сплошная — экспериментальному NjfN. Очевидно, что если числа х,- принадлежат псевдослучай­ ной квазиравномерно распределенной последовательности, то при достаточно боль­ ших N экспериментальная гистограмма (ломаная линия на рис. 4.11, а) приблизится к теоретической прямой 1/т.

Оценка степени приближения, т. е. равномерности последовательности {х*}, может быть проведена с использованием критериев согласия. На практике обычно принимается т —20-Т-50, jV=(102-f-103)/w.

Суть проверки равномерности по косвенным признакам сводится к следующему. Генерируемая последовательность чисел {х,*} разбивается на две последовательности:

Рис. 4.11. Проверка равномерности последователь­ ности

124

^4» ^6» "■) -*2Ь

Затем проводится следующий эксперимент. Если выполняется условие

то фиксируется наступление некоторого события и в счетчик событий добавляется единица. После NJ2 опытов, когда генерировано N число, в счетчике будет некоторое число k ^ N /2 . ___

Геометрически условие (4.13) означает, что точка (х2,_ь х2д, i - 1, N, находится внутри четверти круга радиусом г= 1, что иллюстрируется рис. 4.11, б. В общем случае точка (*2f-i, х2() всегда попадает внутрь единичного квадрата. Тогда те­ оретически вероятность попадания этой точки в четверть круга

P i ~ $1 /4 хруга/^хвадржта= (ЯГ2/4)/(1 • 1)= я /4 .

Если числа последовательности {х,} равномерны, то в силу закона больших чисел теории вероятностей при больших N относительная частота 2k(N-*n(4.

Проверка стохасттности последовательностей псевдослучайных чисел {х/} на­ иболее часто проводится методами комбинаций и серий [7,11,25]. Сущность метода комбинаций сводится к определению закона распределения длин участков между единицами (нулями) или закона распределения (появления) числа единиц (нулей) в л- разрядном двоичном числе Х\. На практике длину последовательности N берут достаточно большой и проверяют все п разрядов или только I старших разрядов числа Xj.

Теоретически закон появленияj единиц в / разрядах двоичного числа Х{ описыва­ ется исходя из независимости отдельных разрядов биномиальным законом рас­ пределения:

р <], 0 - <?/(!)[1 - Р(1)]''У= С^ '(1),

где Р(/, I) — вероятность появленияj единиц в / разрядах числа Х& р (1)=р(0)=0,5 — вероятность появления единицы (нуля) в любом разряде числа Xf, С^=/!/[/!/(/—у)!].

Тогда при фиксированной длине выборки N теоретически ожидаемое число появления случайных чисел Х\ с j ед и н и ц ам и в проверяемых / разрядах будет равно

« ,= № > '( 1).

После нахождения теоретических и экспериментальных вероятностей P{j, I) или чисел ftj при различных значениях /< л гипотеза о стохасгичности проверяется

с использованием критериев согласия [7, И , 18, 21].

При анализе стохастичности последовательности чисел {х,} методом серий по­ следовательность разбивается на элементы первого и второго рода (а и Ъ), т. е.

если Х[<р\

Xf=

в противном случае,

гдеО <р<.

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

После разбиения последовательности {х,} на серии первого и второго рода будем иметь, например, последовательность вида

...aabbbbaaabaaaabbbab... .

Так как случайные числа а и б в данной последовательности независимы и при­ надлежат последовательности (х,}3равномерно распределенной на интервале (0, 1),

125

то теоретическая вероятность появления серии длиной j в последовательности дли­ ной / в N опытах (под опытом здесь понимается генерация числа Xj и проверка условия х i<p) определится формулой Бернулли:

P(j, *)=&{Р (1~Р) J,j = 0 ,1, /=1» п.

В случае экспериментальной проверки оцениваются частоты появления серий длинойj. В результате получаются теоретическая и экспериментальная зависимости P(j\ I), сходимость которых проверяется по известным критериям согласия, причем проверку целесообразно проводить при различных значенияхр, 0 < р < \ и /.

Проверка независимости элементов последовательности псевдослучайных квази­ равномерно распределенных чисел проводится на основе вычисления корреляцион­ ного момента [4].

Случайные величины £ и ц называются независимыми, если закон распределения каждой из них не зависит от того, какое значение приняла другая. Таким образом, независимость элементов последовательности {*,} может быть проверена путем введения в рассмотрение последовательности {у/}= {х/+т}, где т — величина сдвига последовательностей.

В общем случае корреляционный момент дискретных случайных величин £ и ц с возможными значениями д:,-и yj определяется по формуле

i J

где ри — вероятность того, что (£, if) примет значение (х,-, yj).

Корреляционный момент характеризует рассеивание случайных величин f и г\ и их зависимость. Если случайные числа независимы, то К^п= 0. Коэффициент корреляции

Р$ц —KfrjKPxGy)*

где ох—су — средние квадратические отклонения величин f и ij.

При проведении оценок коэффициента корреляции на ЭВМ удобно для вычисле­ ния использовать следующее выражение:

y/D[x&B[Xi+j)

где

2

D[xi+J =

При вычислениях сначала рационально определить суммы:

Xj, ^ Xj+n

126

ется минимальный номер z=z3, при котором истинно событие Р" {хi=xP+i}, и вычисляется длина отрезка апериодичности L —i3+P (рис. 4.12, а). Если Р' оказывается истинным лишь для z= V, то L> V (рис. 4.12, б).

В некоторых случаях достаточно громоздкий эксперимент по определению длин периода и отрезка апериодичности можно заме­ нить аналитическим расчетом, как это показано в следующем при­ мере.

Пример 4.5. Необходимо показать, что в последовательности чисел {*/}, описыва­ емой уравнением

ЛГ|

X i+1s АЛГ/(modМ), */= — ,

при простом модуле М можно так выбрать коэффициент А, что при любом Х0, взаимно простом с М, длина отрезка апериодичности, совпадающая в этом случае с длиной периода Р, будет L = P = M —\. Иначе говоря, надо найти, при каких условиях равенство

(здесь существенно, что Х 0 взаимно просто с М).

По условию требуется, что наименьший показатель степени s= P M(A), удовлет­ воряющий (4.15) и называемый показателем А по модулю М, был равен М —1. Для любого простого модуля М существует ф{М—1) значений А (первообразных корней), удовлетворяющих уравнению (4.15) при Ры (А)= (р (М), где ф(М) — функция Эйлера, определяемая как число натуральных чисел т ^ М , взаимно простых с М. Для простого модуля М имеем ф{М )=М —\.

Таким образом, доказано существование многих А, при которых повторение элементов последовательности {xj наступит на (АГ— 1)-м числе xM_i, т. е. L = P = M —1, что и требовалось доказать.

Для алгоритмов получения последовательностей чисел {JC/} об­ щего вида (4.10) экспериментальная проверка является сложной (из-за наличия больших Р и L), а расчетные соотношения в явном виде не получены. Поэтому в таких случаях целесообразно провести теоретическую оценку длины отрезка апериодичности последовате­ льности L. Для этого воспользуемся элементарной вероятностной моделью, рассмотренной в следующем примере [4, 36, 37].

Пример 4.6. Пусть имеется конечное множество, содержащее N различных чисел. Проведем последовательность независимых опытов, в каждом из которых из этого множества извлекается и записывается одно число. Вероятность извлечения любого числа в каждом из опытов равна X/N, так как выборка чисел проводится с возвратом. Обозначим через L случайную величину — номер опыта, в котором впервые будет снова извлечено уже записанное ранее число. Можно доказать, что в данной вероят­ ностной модели для любого х> 0 имеем

128

lim P { L / ^ N < X } = \ — Q

JV-oo

Так как математическое ожидание случайной величины с такой функцией рас­

пределения равно yjnjl, то при N-+co получим M[L]=y/nNJ2. Такая оценка длины отрезка апериодичности «груба», но полезна на практике для цредварительного определения L с целью дальнейшего уточнения экспериментальным путем.

Рассмотрим некоторые особенности статистической проверки сгохастичности псевдослучайных последовательностей. Для такой проверки могут быть использованы различные статистические кри­ терии оценки, например критерии Колмогорова, Пирсона и т. д. Но в практике моделирования чаще всего пользуются более простыми приближенными способами проверки [29, 37].

Для проверки равномерности базовой последовательности слу­ чайных чисел x h i= l, N, можно воспользоваться такими оценками:

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

Для проверки аппаратных генераторов случайных чисел можно использовать те же приемы, что и для проверки последователь­ ностей псевдослучайных чисел, полученных программным спосо­ бом. Особенностью такой проверки будет то, что проверяются не те числа, которые потом будут необходимы для моделирования систе­ мы S. Поэтому кроме проверки качества выдаваемых генератором случайных чисел должна еще гарантироваться устойчивая работа генератора на время проведения машинного эксперимента с моде­ лью Мм.

Улучшение качества последовательностей. В силу рассмотренных преимуществ основное применение в практике имитационного мо­ делирования систем находят различные программные способы по­ лучения чисел. Поэтому рассмотрим возможные методы улучшения качества последовательностей псевдослучайных чисел. Одним из наиболее употребительных методов такого улучшения является употребление вместо формул вида (4.9), представляющих собой рекуррентные формулы первого порядка, рекуррентных формул порядка г, т. е.

*+1==Ф (*, *_ 1, •••, * —r+i)»

гораздо больше, чем при г=1. Однако при этом возрастает слож­ ность метода, что приводит к увеличению затрат машинного време­ ни на получение чисел и ограничивает возможности его применения на практике.

Для получения последовательности псевдослучайных чисел с бо­ льшой длиной отрезка апериодичности L можно воспользоваться методом возмущений [29, 37]. В основу этого метода получения последовательности чисел положена формула вида

где функции Ф(м) и Ч? (и) различны.

Вэтом случае в основном используется формула х,+1=Ф(х,),

итолько когда i кратно М, последовательность «возмущается», т. е. реализуется переход к формуле х1+ =Ч?(х,). Целое число М на­

зывается периодом возмущения.

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

4.4. МОДЕЛИРОВАНИЕ СЛУЧАЙНЫХ ВОЗДЕЙСТВИЙ НА СИСТЕМЫ

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

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

130