CAPITOLUL 3

 

3.1.2. Metode de generare a numerelor aleatoare

 

       Principalele clase de metode de gene­rare a numerelor aleatoare sunt: metode manuale, metode fizice, metode de memorizare, metode analitice.

 

       Metodele manuale folosesc diferite dispozitive ca: zaruri, urne cu bilete, rulete etc. Un procedeu cunoscut constă în utilizarea unui cilindru prismatic omogen cu 10 feţe numerotate de la 0 la 9. Se aruncă succesiv printr-un mecanism, aleator acest cilindru marcat, obţinându-se aleator una din feţele acestuia. Prin acest procedeu se obţin de fapt numere pseudoaleatoare deoarece intervine neomogenitatea perfectă a cilindrului, uzura în timp a cilindrului etc. care face ca în timp numerele să nu mai aibă probabilităţi egale de apariţie. Metodele manuale se folosesc rar în simularea numerică datorită vitezei reduse.

 

       Metodele fizice se bazează pe analogii dintre unele procese fizice intrinsec întâmplătoare (procese radioactive, procese electronice generatoare de zgomot alb etc.).

       Se poate genera un şir de numere U₁, U₂, ... , U_n, ..., unde U Î (0,1) după cum un detector de particole radioactive înregistrează într-o perioadă determinată de timp Δt, un număr par sau impar de particole emise de sursă. Se ştie că probabilitatea ca să se detecteze k particole în intervalul de timp Δt este dată de legea Poisson

       Probabilitatea ca să se obţină un număr par de particole este

 

       Aceasta înseamnă că trebuie să minimizăm valoarea e^-2λΔt dacă dorim ca 0 şi 1 să apară cu probabilităţi egale, după un număr dat de paşi, cu o eroare oricât de mică dorim.

       Dacă de exemplu dorim ca abaterea e^-λΔt să fie sub 10^-3 vom avea

e^-2λΔt = 10^-3 de unde rezultă că -2λΔt = -3ln10 şi -3(2,3) Þ λΔt » 3,45

       Deoarece λΔt =  rezultă că trebuie ales intervalul de timp Δt astfel încât să se detecteze în medie cel puţin 4 particole.

Aceste şiruri de numere satisfac în cea mai mare măsură caracterul aleator, dar au dezavantajul că nu sunt reproductibile.

 

       Metode de memorizare. O tabelă cu peste 40000 de numere aleatoare "luate la întâmplare din rapoartele de cens" a fost publicată în 1927 de către L.H.C. Tippett. Apoi au fost concepute şi unele mecanisme care să genereze în mod mecanic numere aleatoare. Prima maşină de acest gen a fost folosită în 1939 de către M. G. Kendall şi B. Babington-Smith pentru realizarea unei tabele de 100000 de numere aleatoare. În anul 1955, RAND Corporation a publicat o tabelă larg întrebuinţată, care cuprindea un milion de cifre obţinute cu ajutorul unui alt dispozitiv special. O faimoasă maşină de generare a numerelor aleatoare (numită ERNIE) a fost folosită pentru extragerea numerelor câştigătoare din cadrul loteriei British Premium Savings Bonds.

       Metodele de memorizare a numerelor aleatoare (tabelele) folosesc de regulă memoria internă sau externă a calculatoare­lor. Oferă avantajul reproductibilităţii.

       Prin combinarea metodelor fizice cu memori­zarea pe discuri sau benzi se obţin rezul­tate bune din punct de vedere al preciziei, dar timpul de rulare al programului creşte datorită duratei relativ mare de acces la fiecare înregistrare din memoria externă.

       Se pot folosi şi alte tabele ca de exemplu: cartea de telefoane, tabelele de logaritmi etc. dar numai după efectuarea unor teste de uniformitate

 

       Procedee analitice, folosesc relaţii de recurenţă. Un şir de numere (U_n)_{nÎ N*} se numeşte şir de numere pseudoaleatoare dacă există o anumită formulă de recurenţă şi un număr natural fix k astfel încât

                 U_n = f(U_n-1 ... , U_n-k), n ³ k + 1

       Şirul (U_n),  1 £ n £ k  astfel obţinut, posedă proprietăţile statistice ale unui şir de valori (alese independent) ale unei variabile aleatoare uniform repartizate sau supusă unei anumite legi de repartiţie.

       Există diverse procedee analitice de obţinere a numerelor pseudoaleatoare uniform repartizate, procedee care au fost algoritmizate şi programate pentru calculatoarele electronice existente.

       Un neajuns al acestor algoritmi este că, în aplicarea concretă pe calculator, conduc după un anumit număr de generări la şiruri periodice. Aceasta înseamnă că într-un şir foarte mare de numere pseudoaleatoare există h numere U₁, U₂, .... , U_h cu proprietatea U_i ¹ U_j pentru orice i şi j aparţinând mulţimii {1,2,...,h}. În continuare, însă U_h+1 = U₁, U_h+2 = U₂ etc.

Dacă se repetă la un moment dat s numere (s fiind numărul de valori iniţiale), evident că se repetă întregul subşir de h numere.

       Numărul h reprezintă lungimea intervalului de aperiodicitate. Lungimea perioadei este h-s. Pentru a nu se obţine rezultate eronate este necesar ca şirul generat de numere aleatoare necesar simulării, să nu depăşească lungimea perioadei (h-s), adică: N < h-s unde N reprezintă numărul maxim de cicluri necesare efectuării simulării, în ipoteza că în fiecare ciclu se foloseşte în calcule un singur număr aleator. Dacă această restricţie nu este îndeplinită, este necesar să se considere un alt sistem de s numere iniţiale (neincluse în aceeaşi secvenţă în şirul generat) şi să se genereze un alt şir de numere pseudoaleatoare. Operaţia se repetă până când se obţin cele N numere necesare.

       Una din cele mai cunoscute metode de generare a numerelor pseudoaleatoare este metoda propusă în anul 1946 de John von Neumann "a mijlocului pătratului unui număr" care ţine seama de următoarele considerente:

       Să presupunem că folosim o reprezentare în baza b a numerelor întregi cu care lucrăm (de obicei baza este 2 sau 10) şi că toate aceste numere au 2a cifre (a = 1,2...)(cel puţin după completarea adecvată cu cifre zero în faţa numărului).

       Fiind dat numărul pseudoaleator întreg U_n, următorul număr pseudoaleator U_n+1 se defineşte după von Neumann ca fiind format din cifrele părţii de mijloc a pătratului lui U_n. Ridicând la pătrat pe U_n se obţine un număr cu 4a cifre, luându-se cele 2a cifre de la mijlocul şirului de 4a cifre ale lui se obţine numărul U_n+1.    

(3.1)

unde [u] partea întreagă a lui u.

       Metoda cere numai o singură valoare iniţială U₀ şi are la bază o relaţie simpla de calcul (3.1). Se recomandă ca 2a = 8 şi că la cel puţin două cifre ale numărului să fie diferite de zero şi anume prima cifră obligatoriu ¹ 0 şi cel puţin una pe la mijlocul numărului.

       Ea este însă o sursă slabă de numere pseudoaleatoare deoarece prin acest procedeu de calcul anumite numere se repro­duc. De exemplu, când a = 2, b = 10 numărul 3792 se reproduce 3792² = 14379264

       Condiţiile de repetabilitate sunt stabilite prin studierea ecuaţiei diofantice

       x² - b^ax = r + b^3a · k

unde ţinând seama de condiţiile de repetabilitate avem:

        Mai mulţi cercetători au studiat această metodă la începutul anilor '50. Lucrând cu numere de patru cifre în loc de 10, G. E. Forsythe a încercat 16 valori de start diferite şi a găsit că 12 dintre ele conduc la şiruri care au ciclul 6100, 2100, 4100, 8100, 6100, .... , în timp ce două dintre ele au degenerat luând valoarea zero. N. Metropolis folosind sistemul de numeraţie binar a arătat că atunci când sunt utilizate numere formate din 20 de biţi există 13 cicluri distincte în care şirul degenerează, dintre care cel mai lung are perioada de lungime 142. Lucrând cu numere de 38 de biţi, Metropolis a obţinut un şir de aproximativ 750000 de numere înaintea apariţiei procesului de degenerare, iar cei 750000 x 38 biţi rezultaţi au îndeplinit testele statistice asupra caracterului aleator. Aceasta arată că metoda "mijlocului pătratului unui număr întreg" poate da rezultate utilizabile, dar este destul de riscant să i se acorde prea multă încredere înainte de a o testa minuţios.

       Knuth [25] a demonstrat că a construi o sursă fiabilă de numere aleatoare nu este o întreprindere uşoară. Numerele aleatoare nu se pot produce cu ajutorul calculatoarelor electronice prin metode alese la întâmplare, ci trebuie ca metodele de generare a acesto­ra să se bazeze pe teorii riguroase.

       Printre metodele recurente de generare a numerelor aleatoare, cele care au fost studiate riguros din punct de vedere teore­tic şi au condus practic la rezultate bune, sunt metodele congruenţiale. Acestea au fost iniţiate de D.H. Lehmer în anul 1949. Aceste metode utilizează teoria claselor de restu­ri şi sunt cele mai răspândite.

       Datorită faptului că într-un calculator un număr real poate fi reprezentat numai cu un anumit număr de zecimale, vom genera de fapt numere întregi x_n cuprinse între 0 şi un număr oarecare dat m şi apoi cu ajutorul relaţiei:

vom obţine numere aleatoare cuprinse între zero şi unu. De obicei m este dimensiunea cuvântului de calculator (adică numărul de valori distincte ce pot fi memorate într-un cuvânt de calculator) şi prin urmare, x_n poate fi considerat drept conţinutul întreg al unui cuvânt de calculator cu punctul bazei poziţionat în dreapta iar u_n poate fi considerat conţinutul aceluiaşi cuvânt cu punctul bazei poziţionat în stânga.

       Naylor [35] (1966) clasifică metodele con­gruenţiale în: metode congruenţiale aditive, multiplicative şi mixte.

 

       Metodele congruenţiale aditive

       Se dau r numere iniţiale: x₁, x₂, ... , x_r  şi se generează numere întregi pseudoaleatoare prin formula recursivă:

x_i º (x_i-1 + x_i-r) (mod m), i Î {r+1, r+2, ..}, m Î N^*, unde m este o constantă întreagă, număr prim foarte mare, sau

       x_r+1 º (x₁ + x₂ + ... + x_r) (mod m)

       x_r+2 º (x₂ + x₃ + ... + x_r+1) (mod m)

       x_r+3 º (x₃ + x₄ + ... + x_r+1 + x_r+2) ( mod m)

       :

       :

       :

 

       În general această metodă dă rezultate slabe.

 

       Metode congruenţiale multiplicative, au la bază sistemul de numere magice : [x₀, a, 0, m] unde:

       m,      - modulul;                        m > 0

       a,       - multiplicatorul;               0 £ a < m

       x₀,      - termenul iniţial;    0 £ x₀ < m

 

       Generarea numerelor pseudoaleatoare consecutive se face după relaţia

       x_i+1 º a x_i (mod m), i Î {2,3,..}

 

       Metode congruenţiale mixte

       Aceasta metodă are la bază sistemul de numere [x₀, a, c, m] unde x₀, a şi m sunt mărimile definite anterior iar c o constantă întreagă. Nume­rele generate folosesc clasele de resturi modulo m după relaţia:

                                                  x_i+1 = (ax_i + c) (mod m), i Î {2, 3, ... , 4}

       Şirul obţinut cu ajutorul acestei relaţii se numeşte şir congruenţial liniar.

       Dacă m=10 şi x₀ = a = c = 7 şirul obţinut va fi 7, 6, 9, 0, 7, 6, 9, 0, ....

       Aşa cum se observă, sinul generat nu este "aleator" pentru orice alegere a lui m, a, c şi x₀, existând principii teoretice şi practice riguroase de alegere corespunzătoare a celor patru numere [25]. Exemplul ilustrează faptul că, întotdeauna, şirurile congruenţiale liniare vor intra într-o buclă.    Această proprietate este comună tuturor şirurilor care sunt de forma x_n+1=f(x_n). Ciclul care se repetă se numeşte perioadă. Un şir suficient de aleator va avea întotdeauna o perioadă relativ mare.

       O generalizare a relaţiei de recurenţă a metodei congruenţiale mixte este:

                                                       x_n+k = (a^kx_n + (a^k - 1).c/b)(mod m)

k ³ 0, n ³ 0, b = a-1.

       Scopul fiecărui generator de numere aleatoare este de a obţine şiruri de numere a căror perioadă să aibă lungime maximă. Următoarea teoremă dă o caracterizare riguroasă a situaţiilor în care poate fi realizată perioada de lungime maximă.

 

       Teorema A

       Şirul congruenţial liniar definit de [x₀, a, c, m] are perioada de lungime maximă dacă şi numai dacă:

       i)  c şi m sunt două numere întregi prime între ele;

       ii) b = a-1 este un multiplu de p, pentru orice număr prim p care divide pe m;

       iii) b este multiplu de 4 dacă m este multiplu de 4.

       Demonstraţia teoremei este dată în lucrarea "Tratat de programare a calculatoarelor" de Donald E. Knuth, Editura tehnică, Bucureşti 1983.

       În cazul generatorilor multiplicativi (c = 0) procesul de generare al numerelor aleatoare este mai rapid, lungimea maximă a perioadei numerelor generate este relativ mică.

       Realizarea unei perioade suficient de mari, presupune îndeplinirea anumitor condiţii de către multiplicatorul a.

       Fie λ(m) ordinul unui element primitiv, adică maxim posibil modulo m.

Dacă a şi m sunt două numere prime între ele, atunci cel mai mic număr întreg λ pentru care a^λ  º 1 (modulo m) este numit, prin convenţie ordinul lui a modulo m. Orice astfel de număr a care are ordinul modulo m maxim posibil este numit element primitiv modulo m.

Din teorema lui Euler rezultă că λ(p^e) este un divizor al numărului p^e-1•(p-1). În acest caz λ(2) = 1, λ(4) = 2, λ(2^e) = 2^e-2 dacă e ³ 3; λ(p^e) = p^e-1(p-1), dacă p > 2,     λ(p₁^e1 ... p_i^ei) = c.m.m.m.c. (λ(p₁^e1), ... , λ(p_i^ei))

 

       Teorema B [R.D. Carmichael].

       Lungimea maximă posibilă a perioadei pentru cazul în care c = 0, este λ(m) definit ca mai sus. Perioada de lungime maximă este realizată dacă:

       i) x₀ şi m sunt prime între ele;

       ii) a este element primitiv modulo m.

       Dacă m este număr prim, atunci putem obţine o perioadă de lungime m-1, această valoare este doar cu o unitate mai mică decât lungimea maximă a perioadei şi, prin urmare, o astfel de perioadă este potrivită oricărui scop practic.

       Modul de construire a unor generatori multiplicativi congruenţiali este prezentat de Văuva Ion  în "Modele de simulare cu calculatorul", Editura tehnică Bucureşti, 1977 şi Knuth E. Donald în "Tratat de programarea calculatoarelor", Editura tehnică, Bucureşti, 1983. Dintre aceştia, un generator perfor­mant care satisface condiţiile necesare este generatorul RND de forma [x₀,16807,0,2³¹-1] modulul m = 2³¹-1 fiind cel mai apropiat număr prim de cuvântul calculatoarelor FELIX C - 256 şi IBM - 360. Acest generator este programat în limbajul ASSIRIS şi este accesibil pentru calculatoarele Felix C - 256. Subrutina RND(IX,YFL) generează o valoare de selecţie YFL uniformă pe (0,1), utilizând ca număr de pornire ("sămânţa") pe IX, 0<IX<2³¹-1 iar apelarea ei este făcută din orice program FORTRAN prin CALL RND (IX,YFL).

       O altă subrutină sau generator de numere aleatoare este RANDU de forma [X₀,1230703125,0,2³¹-1] iar apelarea din programe FORTRAN este realizată prin CALL RANDU (IX, IY,YFL) în care:

 

       -        IX este numărul de pornire (Număr întreg impar cu mai puţin de 9 cifre); numărul IX trebuie definit (recursiv sau nu) înainte de orice apelare a subrutinei;

       -        IY număr întreg obţinut prin RANDU şi are o valoare cuprinsă între 1 şi 2³¹-1. El poate fi utilizat la o nouă apelare a subrutinei când IX = IY din apelarea precedentă;

       -        YFL număr aleator în intervalul (0,1) obţinut de subrutină.

 

       În continuare vom da secvenţa instrucţiunilor FORTRAN pentru generarea a 10 numere aleatoare uniform distribuite în intervalul (0,1) şi instrucţiunile subrutinei RANDU.

       IX = 12341

       DO  1  I=1,10

       CALL RANDU  (IX,IY,YFL)

       IX=IY

       WRITE (3, 200)  YFL

200  FORMAT (10X, F10.8)

  1  CONTINUE

       STOP

       END

       SUBROUTINE RANDU (IX,IY,YFL)

       IY = IX * 1230703125

       IF(IY) 10,20,20

 10  IY = IY + 2147483647 + 1

 20  YFL = IY * 0.4656613E - 9

       RETURN*

       END

Fiecare limbaj de simulare şi aproape fiecare limbaj de programare îşi are propriul său generator de numere aleatoare uniform distribuite.

       În Anexa 1 se prezintă posibilităţile oferite de limbajul C+, pentru generarea numerelor aleatoare uniform distribuite.

       O altă posibilitate de a obţine şiruri pseudoaleatoare uniform distribuite în intervalul [0,1] este oferită de metodele comparative. În acest caz se presupune că dispunem de n numere pseudoaleatoare repartizate uniform în intervalul [0,1] U₁, U₂, ... , U_n. În acest caz numărul U_n+p se obţine după următorul algoritm:

       Pasul 1.

       Se citeşte şirul iniţial: U₁, U₂, ... , U_n (n ³ 100)

       Pasul 2.

       Dacă U_n+p £ U_p+1 atunci U_n+p+1 = U_n+p + (1 - U_n+p)U_p+2

       În caz contrar se trece la pasul 3.

       Pasul 3.

                                                                 U_n+p+1 = U_n+p • U_p+3

       Pasul 4.

       Se reia algoritmul de la pasul 2 până când se ajunge la lungimea dorită a şirului de numere pseudoaleatoare.

       Pentru a îmbunătăţi calitatea numerelor generate, se pot utiliza metodele combina­ţionale de generare a numerelor pseudoalea­toare uniform repartizate.

       De exemplu, un algoritm cu doi genera­tori (G1 şi G2) are următorii paşi:

       Pasul 0.

       Se iniţializează programul cu p numere uniform distribuite în intervalul [0,1]: U₁, U₂, ... , U_p.

       Pasul 1.

       Cu ajutorul lui G1 se generează un număr întreg, pseudoaleatoar, uniform repartizate pe intervalul [1,p].

       Numărul generat n este folosit ca indice pentru a extrage din p numere pseudo­aleatoare uniform repartizate, un număr, care se înregistrează în secvenţe de numere utilizate în modelul de simulare.

       Pasul 2.

       Al doilea generator G2, bazat pe metoda multiplicativă, generează un nou număr uniform repartizat pe [0,1] care înlocuieşte numărul extras la pasul 1, astfel ca să se completeze şirul celor p numere pseudoaleatoare uniform repartizate.

       Pasul 3.

       Se repetă pasul 1 şi 2 ori de câte ori este necesar astfel încât să se obţină lungimea dorită a şirului de numere necesare simulării.