Loading [MathJax]/extensions/tex2jax.js
Logo ROOT   6.12/07
Reference Guide
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
quasirandom.C File Reference

Detailed Description

View in nbviewer Open in SWAN Example of generating quasi-random numbers

pict1_quasirandom.C.png
Processing /mnt/build/workspace/root-makedoc-v612/rootspi/rdoc/src/v6-12-00-patches/tutorials/math/quasirandom.C...
Time for gRandom Real time 0:00:00, CP time 0.000
Time for Sobol Real time 0:00:00, CP time 0.000
Time for Niederreiter Real time 0:00:00, CP time 0.000
number of empty bins for pseudo-random = 31139
number of empty bins for sobol = 30512
number of empty bins for niederreiter-base-2 = 30512
(int) 0
#include "Math/Random.h"
#include "TH2.h"
#include "TCanvas.h"
#include "TStopwatch.h"
#include <iostream>
using namespace ROOT::Math;
int quasirandom(int n = 10000, int skip = 0) {
TH2D * h0 = new TH2D("h0","Pseudo-random Sequence",200,0,1,200,0,1);
TH2D * h1 = new TH2D("h1","Sobol Sequence",200,0,1,200,0,1);
TH2D * h2 = new TH2D("h2","Niederrer Sequence",200,0,1,200,0,1);
// quasi random numbers need to be created giving the dimension of the sequence
// in this case we generate 2-d sequence
// generate n random points
double x[2];
TStopwatch w; w.Start();
for (int i = 0; i < n; ++i) {
r0.RndmArray(2,x);
h0->Fill(x[0],x[1]);
}
std::cout << "Time for gRandom ";
w.Print();
w.Start();
if( skip>0) r1.Skip(skip);
for (int i = 0; i < n; ++i) {
r1.Next(x);
h1->Fill(x[0],x[1]);
}
std::cout << "Time for Sobol ";
w.Print();
w.Start();
if( skip>0) r2.Skip(skip);
for (int i = 0; i < n; ++i) {
r2.Next(x);
h2->Fill(x[0],x[1]);
}
std::cout << "Time for Niederreiter ";
w.Print();
TCanvas * c1 = new TCanvas("c1","Random sequence",600,1200);
c1->Divide(1,3);
c1->cd(1);
h0->Draw("COLZ");
c1->cd(2);
// check uniformity
h1->Draw("COLZ");
c1->cd(3);
h2->Draw("COLZ");
gPad->Update();
// test number of empty bins
int nzerobins0 = 0;
int nzerobins1 = 0;
int nzerobins2 = 0;
for (int i = 1; i <= h1->GetNbinsX(); ++i) {
for (int j = 1; j <= h1->GetNbinsY(); ++j) {
if (h0->GetBinContent(i,j) == 0 ) nzerobins0++;
if (h1->GetBinContent(i,j) == 0 ) nzerobins1++;
if (h2->GetBinContent(i,j) == 0 ) nzerobins2++;
}
}
std::cout << "number of empty bins for pseudo-random = " << nzerobins0 << std::endl;
std::cout << "number of empty bins for " << r1.Name() << "\t= " << nzerobins1 << std::endl;
std::cout << "number of empty bins for " << r2.Name() << "\t= " << nzerobins2 << std::endl;
int iret = 0;
if (nzerobins1 >= nzerobins0 ) iret += 1;
if (nzerobins2 >= nzerobins0 ) iret += 2;
return iret;
}
Author
Lorenzo Moneta

Definition in file quasirandom.C.