|
Для моделирования ста |
|
|
|
|
|||||||
ционарных процессов весь |
|
|
|
|
||||||||
ма |
эффективными |
оказа |
|
|
|
|
||||||
лись |
методы, |
оспованные |
|
|
|
|
||||||
на |
применении |
быстрого |
|
|
|
|
||||||
преобразования |
|
Фурье |
|
|
|
|
||||||
(БПФ) |
1451. |
Для |
|
модели |
|
|
|
|
||||
рования |
гауссовского ста |
|
|
|
|
|||||||
ционарного |
|
случайного |
|
|
|
|
||||||
процесса (ГССП) по задан |
|
|
|
|
||||||||
ной |
спектральной |
плот |
Рис. 42. Блок-схема алгоритма моде |
|||||||||
ности |
W (/) |
в |
работе [86] |
лирования амплитуды a h к -го цикла |
||||||||
предложен |
алгоритм, |
ос |
с |
помощью равномерно-распределен |
||||||||
нованный на модификации |
ной псевдослучайной величины а. |
|||||||||||
известного |
разложения |
в |
|
|
|
At ординаты |
||||||
ряд Фурье |
[26]. |
Равноотстоящие по |
времени |
|||||||||
процесса хтвычисляются по формуле |
|
|
||||||||||
Хт = |
£ . F*-4*exp(/Jj£- mfcV |
т — 0, |
1, |
N — i, (2.79) |
||||||||
|
|
/г=0 |
|
|
|
|
|
|
|
|
|
|
где N — число отсчетов в моделируемой реализации случай ного процесса; f = —1; F/, — комплексная гауссовская
случайная величина:
Vh = V'N- h, |
|
E(Vh V9l) = 8lh] |
(2.80) |
||
V* — комплексно-сопряженная |
V\ |
Е (•) — символ ос |
|||
реднения; bih — 1 при I = |
к |
и bih = |
0 при I ^ |
к. Коэффи |
|
циенты Ли вычисляются следующим образом: |
|
||||
A„ = V W (Afk); |
* = 0, |
1, . . . , 4 — |
1; |
||
-4а = A N —^ |
к = |
- у - » • • • » N — 1; |
(2.81) |
||
где |
|
1 |
|
|
|
Д/ = |
' |
|
(2.82) |
||
N A t |
|
||||
Формула (2.82) связывает шаг по частоте А/, шаг по вре мени АЬи число точек отсчета N. Для корреляционной функ
ции случайной последовательности хт, используя выражения (2.79) — (2.81), докажем такое равенство:
R (n\t) = Е (хт+пхт) = “S W (Д/Л) exp (/ |
nkj . (2.83) |
Действительно,
Е (^m+Fv^ro) —■
= |
Е |
S |
VkA hV*A*i exp |/ |
(— in k + m i + |
ш)jj = |
|
N |
- l JV—1 |
г |
|
1 |
= |
S |
S ^(VfeF[)^/,i4iexp |
m/c + mi + |
ш) = |
|
|
A=0 i=0 |
L " |
|
J |
|
= 2 A l exp f j -4P- ra/c |
|
h=0 |
N |
|
|
Из выражений (2.79) и (2.83) следует, что массивы чисел Хт и R (nAt) получаются из массивов VuAh и W (А/к) соот
ветственно с помощью обратного дискретного преобразования Фурье (ОДПФ) [45]. По свойствам преобразования Фурье с учетом равенств (2.80) и (2.81) следует, что хт и R (nAt) —
действительные числа. Изложенное выше позволяет считать, что Хт при достаточно большом N представляют собой дис
кретные отсчеты ГССП с нулевым средним значением и спект ральной плотностью W (/). Среднеквадратическое значение
процесса
а , = V Е (х2т) = К 5 Щ = ( У |
Л,:)и |
(2.84) |
'А=0 |
/ |
|
Из выражений (2.79), (2.81) и (2.82) следует, что наивысшая частота при воспроизведении спектра процесса
AfN _ 1
(2.85)
2 2Д*
Это значит, что число моделируемых точек на период гар моники с максимальной частотой AfN/2 равно 2. Для опре
деления экстремумов процесса необходимо, чтобы число точек на период было не менее 8. Поэтому, если спектральная плот ность процесса W (/) задана на интервале 0 ^ ^ fm (fm — максимальная частота спектра), то А/ и At нужно принимать
такими:
At = -gjr-; д 1 = |
(2.86) |
При этом, применяя формулу (2.81), при А/к > |
fm полагают, |
что W (Afk) = 0.
Последовательность хь получается из последовательности VkAh с помощью ОДПФ. Если N — 2т (т — натуральное
число), то для выполнения ОДПФ можно применить алго ритм БПФ [45]. Алгоритм БПФ позволяет существенно, при мерно в Nl\og2N раз сокращать время счета по сравнению с прямыми вычислениями по формуле (2.79), что для N — 4096
составляет больше двух порядков. При этом результат БПФ формируется в том же массиве, в котором были заданы коэф
фициенты VhAk. Поскольку УкАц — комплексные величины, то массив {VhAk} — комплексный, а результат преобразо
вания Фурье — массив {#m} действительный. Поэтому числа хт содержатся в действительной части исходного массива. Практически без изменения объема вычислений можно видо изменить алгоритм так, чтобы моделировать сразу две реали зации случайного процесса, одпа из которых будет форми роваться в действительной части исходного массива, а вто рая — в мнимой части. Для этого модифицируем (2.79) так:
|
|
N— 1 |
/ |
\ |
|
|
Zm= |
Хт+ |
)Ут= 2 |
У*Лнexpf/ -Щ- тк\ + |
|
||
+ Д jTkA he |
x p |
mkj = |
Д |
+ Л ) Ahexp |
mkj |
|
(2.87)
Чтобы хти утбыли отсчетами требуемой реализации ГССП, необходимо, чтобы каждый из массивов гауссовских случай ных величин (F/,} и {Т^} удовлетворял соотношению (2.80).
Для реализации алгоритма (2.87) необходимо предложить способ формирования массива случайных комплексных ве
личин Uk = |
Vi, + ;Т»,. Покажем, |
что для этого достаточно |
||||
вычислять величины Uk так: |
|
|
|
|||
|
Uh = uh -\-jwh\ |
к = 0, |
1, . . . , |
N — 1, |
(2.88) |
|
где ик и |
wu — действительные |
взаимно-независимые гаус- |
||||
совские |
случайные величины, |
(щ) = |
(wk) = 0, |
(ик) — |
||
= {и>п) = |
1. |
Случайные |
числа |
ику wk |
моделируются на |
|
ЭВМ с помощью стандартных подпрограмм. Рассмотрим мас сивы случайных величин:
|
|
|
|
|
|
|
wh ~ WN - h |
(2.89) |
|
|
|
|
|
|
|
|
|
|
|
|
т |
|
“>h + |
u*N—k |
, , |
aN—h ~ |
“ft |
(2.90) |
|
|
i ft = |
------ о-------- г 7 |
|
|
|||||
|
|
|
|
uh+ |
UJV—ft . |
Wu —W |
|
||
Vk + jTk |
= |
. “'ft |
N - h |
||||||
|
|
+ )■ |
|
+ |
|||||
, |
wN - h |
+ wN - ( N - k ) |
uJV-ft “ |
uft |
_ n |
||||
- t ----------- 2------------------------ 2---------- ° ht |
|||||||||
Покажем, |
что |
Vk и |
Тк |
удовлетворяют |
условиям (2.80): |
||||
|
aN - h |
+ uN - ( N - h ) 2 |
ID* |
wN —(N—k) |
|||||
|
’N - h ~ |
||||||||
|
|
|
2 |
|
|
— ] |
|
|
|
|
uk + |
UN—ft |
+ |
wh |
wN—h |
= |
Vhi |
||
|
|
|
|
|
/ |
2 |
|||
Е (УhVh) = |
[Е (Uft 4- Mjv—ft)2 + E (Wh— H>jv_ft)2] = |
— ~£~\E (щ) -f- E (Wh) + E (u%—h) -f- E (u?^_;t)] = 1;
|
s {VkVi-k) = E |
|
j |
x |
||
|
|
x |
^ «*iV-fc + *ft |
. wN - h ~ Wh у = |
|
|
1 |
|
|
[(uk + MW_ft) { w N - k — |
|
||
= T |
E [(“ ft + “ JV—ft)2] ------J - E |
Wft)] |
||||
|
|
+ |
~J~ E l(u>ft — WN—k) (UN—k + Wft)] + |
|
|
|
|
|
1 |
|
|
|
|
|
У |
— |
Е [ { w k — I V N - k ) |
{ w N - h — Wft)] = |
0. |
|
Аналогично |
доказывается, что |
E (FfcF*) = 0 |
при |
l Ф k н |
||
l Ф N |
— к. |
Доказательство |
справедливости |
соотношений |
||
(2.80) |
для Тл идентично доказательству для |
F*. Докажем, |
||||
что Vh ж Th независимы. |
|
|
|
|||
Случай |
I: к Ф 1 7 1 ф N — к. |
Доказательство, очевидно, |
следует из |
того, что в формулах |
(2.89) и (2.90) фигурируют |
независимые случайные величины. |
|
|
||||
Случай И: к = |
I. Тогда |
|
|
|
||
Е (VhT'h) = |
Е [(■ “А+ «Л_Л |
, |
. U>h- wN-h |
\ |
||
|
|
|
~2-------+ >--------2-------) |
|||
|
|
+ Ц”4 |
~ Вй )] = |
4 - Я [(«. + |
»*-*) X |
|
X {wk + |
WN—а)] + |
Е [(Ufe + |
UJV—л) (UN—A— Hft)] -j- |
|||
+ |
4 ~ E [ { w h — |
W N - k ) |
{ w h |
+ W N - ft)] — |
|
|
|
l |
|
|
(tiff—ft — Uft)] = |
|
|
-------^ - E [ ( w h — wN- k ) |
|
|||||
= 4 “ \M iuN—h) — E (ul) + E (wl) — E (wjv_ft)] — 0.
Случай III: к = N — L Доказательство проводится ана
логично доказательству случая II.
Таким образом, если вычислять Vh согласно формуле (2.88), то после применения ОДПФ к массиву (UhAh) по
формуле (2.87) получим массив комплексных чисел {zm}, действительная и мнимая части которых представляют собой
отсчеты Двух независимых ГССП с заданной с помощью мас сива {Л/,} спектральной плотностью. Для оптимизации алго ритма и программы моделирования случайной нагрузки не обходимо также учесть следующие два момента. Во-первых, целесообразно моделировать ГССП с единичной дисперсией, поэтому с учетом выражения (2.84) в соотношении (2.87) сле дует использовать нормированную спектральную плотность, т. е. величины 4 “ вычислять так:
Во-вторых, для реализации ОДПФ методом БПФ можно пользоваться программами, которые вычисляют прямое ДПФ методом БПФ. Действительно, формулу (2.87) можно запи сать так:
z [ S ^ e x p f ; - ^ ] ]
[ , ? п # JUft ехР ( — / ~ т г тЛ :) ] |
(2.92) |
|
Поскольку случайные величины 7д [см. формулу (2.90)] и ут [см. формулу (2.87)] являются центрированными, то с ве роятностной точки зрения £/д и С/д, а также Zm и zm нераз
личимы. Поэтому для реализации алгоритма моделирования пригодны без всяких изменений программы, вычисляющие БПФ.
После расчета реализации необходимо выделить экстре мумы нагрузки. Опишем алгоритм выделения экстремума на примере максимумов. При обработке реализации просмат ривались тройки последовательных отсчетов: дг*_хс, дгц-i.
Считали, что максимум достигается в точке it если выполня лись условия Xi—i <С Xi tfi+i. При этом, однако, действи
тельное значение максимума несколько больше, чем по скольку точка дискретного отсчета может быть смещена от точки максимума. Поэтому по значениям Xi—i; xt; Xi+i на
ходили значение соответствующего максимума так: прово дили квадратичную параболу (рис. 43) через точки (t — At; Xi—i); (t; Xi); (t + At; ari+i) и находили наибольшее значение функции х:
х = xi -f- 0,125 (яц.1 — Xi—\)2f(2xi — — ач-i), (2.93)
которое и принимали за максимум.
Расчеты показали, что погрешность формулы (2.93) при нахождении максимума синусоиды но дискретным отсчетам