KADD 2022 Laboratorium 6 EN

From Łukasz Graczykowski

---

Revision as of 11:26, 4 April 2022

Contents 1 Exercise 2 Uwagi 3 Wynik

Exercise

Part one: linear congruent generator of pseudorandom numbers (1 pkt.)

Please write a generator of pseudorandom numbers and save generated numbers into a file.

The generator should be based on the following formula:

x[j+1] = (g*x[j] + c) mod m.

This generator is called LCG - linear congruent generator. Providing the first (beginning) number x[0] defines the whole series, which is periodic. The period depends on the parameters and under certain conditions it reaches maximum value of m. These conditions are:

• c and m do not have joint divisors,
• b = g-1 is a multiply of every prime number p, which is a divisor of a number m,
• b is a multiply of 4 if n is also a multiply of 4.

For simplicity we can uuse c = 0, and in that case we get a multiplicative generator (MLCG).

• Value of g and m should be easy to modify in the program.

The result of the macro execution should be a file with the name name.dat containing a series of generated numbers for a given set of parameters. The macro should be executed three times, resulting in three files: random1.dat, random2.dat, random3.dat, for the following parameters, respectively:

• m=97 i g=23,
• m=32363 i g=157,
• m=147483647 i g=16807.

Część druga: test widmowy (1 pkt.)

Należy przeprowadzić test widmowy aby przetestować jakość generatora. By to zrobić należy narysować na płaszczyźnie punkty o współrzędnych (x[n], x[n+1]). Uzyskany obraz utworzy wzór przypominający widmo generatora - stąd nazwa testu.

Jeśli punkty będą rozłożone równomiernie generator można uznać za dobry. Jeśli zdecydowanie widać pewną okresowość - punkty powtarzają się wielokrotnie - generator nie działa poprawnie. Oczywiście na rozłożenie punktów wpływa jedynie dobór parametrów g i m.

• Do tworzenia wykresów widma poleca się użycie obiektów TH2D.

Wynikiem powinny być trzy wykresy widma.

Część trzecia: generacja liczb losowych oparta na transformacji rozkładu jednorodnego (3 pkt.)

Dowolna funkcja zmiennej losowej jest zmienną losową. Powstaje więc pytanie jaka jest gęstość zmiennej losowej Y jeżeli znana jest gęstość f(x). Zakładamy, że prawdopodobieństwo g(y)dy jest równe f(x)dx, gdzie dx odpowiada wartością dy. Warunek jest spełniony dla odpowiednio małych dx. Wynika stąd, że:

g(y) = dx/dy f(x)

Teraz jeżeli założymy, że gęstość prawdopodobieństwa f(x) wynosi 1 w 0<=x<=1 i f(x) = 0 dla x<= 0 i x>1 to powyższe równanie możemy zapisać w postaci:

g(y)dy = dx = dG(y),

gdzie G(y) jest dystrybuantą zmiennej losowej Y. Co po całkowaniu daje nam

x = G(y) => y = G^-1(x).

Jeśli zmienna losowa X ma rozkład jednostajny na odcinku pomiędzy 0 i 1 oraz jeśli znana jest funkcja odwrotna G^-1(x) to funkcja g(y) opisuje gęstość prawdopodobieństwa rozkładu zmiennej losowej Y.

Używając tej metody należy wygenerować 10000 liczb z rozkładu:

Dla tau = 2:

• Należy wygenerować 10000 liczb z rozkładu 0 do 1 używając generatora z części pierwszej (zapisać do pliku wartości xn makrem z pierwszej części a następnie je wczytać w makrze z części drugiej).
• Analitycznie (na kartce) policzyć dystrybuantę tego rozkładu, a następnie funkcję odwrotną. (1 pkt.)
• Wygenerować rozkład f(x) - wrzucając wygenerowane wartości do histogramu - korzystając z: (1 pkt.)
  • liczb wygenerowanych wcześniej i wczytanych z plików losowe1.dat, losowe2.dat, losowe3.dat,
  • standardowego generatora ROOT'a, np. gRandom->Uniform(1) (obiekt gRandom istnieje domyślnie w uruchomionej instancji ROOT; można oczywiście również stworzyć samodzielnie obiekt TRandom - link).
• Narysować na jednym wykresie histogram (odpowiednio unormowany) oraz funkcję teoretyczną f(x) (obiekt TF1). (1 pkt.)

Uwagi

• Czytamy dokładnie Wykład 4 (link), zwłaszcza slajdy 6 oraz 18-25
• Na samym początku, przed losowaniem, musimy samodzielnie ustawić wartość pierwszej liczby pseudolosowej x0 (tzw. ziarno, "seed"). Jeżeli chcemy, by za każdym razem liczby pseudolosowe były inne, możemy je ustawić z zegara systemowego:

  x0 = time(NULL);

• Parametry histogramów z obrazków poniżej:

 TH1D *hUniform = new TH1D("hUniform","Uniform distribution",100,0,1);
 TH2D *hCorr = new TH2D("hCorr","Correlation",100,xmin,xmax,100,0,1);

• Ilość losowań w części pierwszej:

 const int N = 1000000;

• Wczytywanie danych z pliku:

ifstream ifile;
ifile.open("dane.dat");
double val;
while(ifile>>val)
{
  cout<<val<<endl;
}
ifile.close();

• Zapisywanie danych do pliku:

ofstream ofile;
ofile.open("dane.dat");
for(int i=0;i<N;i++)
  ofile<<val<<endl;
}
ofile.close();

Wynik

Przykładowy rozkład dla parametrów:

• m=97, g=23

• m=2147483647, g=16807

Przykładowy wynik transformacji rozkładu jednorodnego: