| 1 | // Test des operations de tirages aleatoires thread safe/nosafe | 
|---|
| 2 | #include "sopnamsp.h" | 
|---|
| 3 | #include "zthread.h" | 
|---|
| 4 |  | 
|---|
| 5 | #include <stdlib.h> | 
|---|
| 6 | #include <stdio.h> | 
|---|
| 7 | #include <unistd.h> | 
|---|
| 8 | #include <string> | 
|---|
| 9 | #include <iostream> | 
|---|
| 10 |  | 
|---|
| 11 | #include "tvector.h" | 
|---|
| 12 | #include "matharr.h" | 
|---|
| 13 | #include "fioarr.h" | 
|---|
| 14 | #include "tarrinit.h" | 
|---|
| 15 | #include "srandgen.h" | 
|---|
| 16 | #include "nbtri.h" | 
|---|
| 17 | #include "timing.h" | 
|---|
| 18 |  | 
|---|
| 19 | #include "randr48.h" | 
|---|
| 20 |  | 
|---|
| 21 | long Check_Redondances(TVector<r_8>& v); | 
|---|
| 22 |  | 
|---|
| 23 | //--------------------------------------------------------------------- | 
|---|
| 24 | // Classe heritant de ZThread, pour remplir un vecteur avec des aleatoires | 
|---|
| 25 | class MTRndTest : public ZThread { | 
|---|
| 26 | public: | 
|---|
| 27 | MTRndTest(TVector<r_8>& v,bool thsafe=true,int nbuff=1024); | 
|---|
| 28 | virtual void run(); | 
|---|
| 29 | protected: | 
|---|
| 30 | string nom_; int inum_; | 
|---|
| 31 | TVector<r_8> vv_; | 
|---|
| 32 | bool thsafe_; | 
|---|
| 33 | int nbuff_; | 
|---|
| 34 | }; | 
|---|
| 35 |  | 
|---|
| 36 | static int mtrandId = 0;  // Pour donner un identificateur a chaque thread | 
|---|
| 37 |  | 
|---|
| 38 | MTRndTest::MTRndTest(Vector& v,bool thsafe,int nbuff) | 
|---|
| 39 | : vv_(v,true) , thsafe_(thsafe), nbuff_(nbuff) | 
|---|
| 40 | { | 
|---|
| 41 | if(nbuff_<=0) nbuff_ = 1024; | 
|---|
| 42 | char buff[32]; | 
|---|
| 43 | inum_ = mtrandId; | 
|---|
| 44 | sprintf(buff,"MTRndTest-Id=%d",mtrandId); | 
|---|
| 45 | mtrandId++; | 
|---|
| 46 | nom_ = buff; | 
|---|
| 47 | cout<<"Thread MTRndTest("<<nom_<<" ) Created ... with threadsafe="<<thsafe_ | 
|---|
| 48 | <<", nbuff="<<nbuff_<<endl; | 
|---|
| 49 | } | 
|---|
| 50 |  | 
|---|
| 51 |  | 
|---|
| 52 | // Le travail fait dans chaque thread (appelle par start()) | 
|---|
| 53 | void MTRndTest::run() | 
|---|
| 54 | { | 
|---|
| 55 | cout<<"MTRndTest::run() - Nom= "<<nom_<<" vv.Size()= "<<vv_.Size() | 
|---|
| 56 | <<" with threadsafe="<<thsafe_<<endl; | 
|---|
| 57 | if(thsafe_) {  // On appele ThSDR48RandGen Thread-Safe | 
|---|
| 58 | ThSDR48RandGen rgen(nbuff_,true); | 
|---|
| 59 | for(sa_size_t k=0; k<vv_.Size(); k++) vv_(k) = rgen.Flat01(); | 
|---|
| 60 | } else {  // On appele ThSDR48RandGen NO Thread-Safe | 
|---|
| 61 | // ATTENTION: le comportement est different selon qu'un thread | 
|---|
| 62 | // se termine avant ou apres que le suivant commence: | 
|---|
| 63 | // Pour de petits VSZ le run() du thread "n" se termine | 
|---|
| 64 | // avant que le run() du thread "n+1" ne commence: | 
|---|
| 65 | // cela donne un comportement "tread_safe" APPARENT mais | 
|---|
| 66 | // dans la realite drand48() N'EST PAS thread_safe !!!! | 
|---|
| 67 | ThSDR48RandGen rgen(0,false); | 
|---|
| 68 | for(sa_size_t k=0; k<vv_.Size(); k++) vv_(k) = rgen.Flat01(); | 
|---|
| 69 | } | 
|---|
| 70 | cout<<"....MTRndTest::run() - "<<nom_<<" ENDED !!!!!!!!!"<<endl; | 
|---|
| 71 | } | 
|---|
| 72 |  | 
|---|
| 73 | //--------------------------------------------------------------------- | 
|---|
| 74 | int main(int narg, char *arg[]) | 
|---|
| 75 |  | 
|---|
| 76 | { | 
|---|
| 77 | SophyaInit(); | 
|---|
| 78 |  | 
|---|
| 79 | if(narg<3) {cout<<"usage: tmtrnd3 ntir nthreads [nbuff]"<<endl; return -1;} | 
|---|
| 80 | sa_size_t VSZ = atoi(arg[1]); | 
|---|
| 81 | int NTH = atoi(arg[2]); | 
|---|
| 82 | int NBUFF = (narg>3) ? atoi(arg[3]): 64; | 
|---|
| 83 | if(NTH<=0) NTH = 1; | 
|---|
| 84 | cout<<"tmtrnd/arguments, VSZ= "<<VSZ<<" NTH= "<<NTH<<" NBUFF="<<NBUFF<<endl; | 
|---|
| 85 | // Il ne faut pas que les buffers de ThSDR48RandGen soient remplis | 
|---|
| 86 | // au dela de ce qui est necessaire si on veut faire | 
|---|
| 87 | // la comparaison avec le mono thread | 
|---|
| 88 | int lmod = NTH*NBUFF; | 
|---|
| 89 | if(VSZ<lmod) VSZ = lmod; | 
|---|
| 90 | if(VSZ%lmod!=0) VSZ = int(VSZ/(double)lmod+0.5)*lmod; | 
|---|
| 91 | cout<<"tmtrnd/starting, lmod="<<lmod<<" NTH= "<<NTH<<" -> changed VSZ= "<<VSZ<<endl; | 
|---|
| 92 |  | 
|---|
| 93 | int rc = 0; | 
|---|
| 94 | sa_size_t csz = VSZ / NTH; | 
|---|
| 95 | sa_size_t first, last; | 
|---|
| 96 |  | 
|---|
| 97 | try { | 
|---|
| 98 |  | 
|---|
| 99 | //--- | 
|---|
| 100 | //--- Init et set de l'aleatoire | 
|---|
| 101 | //--- | 
|---|
| 102 | ThSDR48RandGen rgen_nothsf(0,false); | 
|---|
| 103 | unsigned short seed_16v[3]; | 
|---|
| 104 | { | 
|---|
| 105 | cout<<"\n>>>>>>> Init de l'aleatoire"<<endl; | 
|---|
| 106 | rgen_nothsf.AutoInit(5); | 
|---|
| 107 | rgen_nothsf.GetSeed(seed_16v); | 
|---|
| 108 | } | 
|---|
| 109 |  | 
|---|
| 110 | InitTim(); | 
|---|
| 111 |  | 
|---|
| 112 | //--- | 
|---|
| 113 | //--- remplissage en mono-thread directement avec drand48 | 
|---|
| 114 | //--- | 
|---|
| 115 | TVector<r_8> DATA0(VSZ); | 
|---|
| 116 | { | 
|---|
| 117 | cout<<"\n>>>>>>> Remplissage en mono-thread avec drand48()"<<endl; | 
|---|
| 118 | rgen_nothsf.SetSeed(seed_16v); | 
|---|
| 119 | for(int i=0;i<VSZ;i++) DATA0(i) = drand48(); | 
|---|
| 120 | PrtTim(" "); | 
|---|
| 121 | } | 
|---|
| 122 |  | 
|---|
| 123 | //--- | 
|---|
| 124 | //--- Remplissage par threads avec ThSDR48RandGen thread safe | 
|---|
| 125 | //--- | 
|---|
| 126 | TVector<r_8> DATA1(VSZ); | 
|---|
| 127 | { | 
|---|
| 128 | cout<<"\n>>>>>>> Remplissage par threads avec ThSDR48RandGen thread safe"<<endl; | 
|---|
| 129 | rgen_nothsf.SetSeed(seed_16v); | 
|---|
| 130 | vector<MTRndTest *> vth1; | 
|---|
| 131 | cout<<"...tmtrnd/creating threads "<<endl; | 
|---|
| 132 | for(int kt=0; kt<NTH; kt++) { | 
|---|
| 133 | first = kt*csz; | 
|---|
| 134 | last = (kt == NTH-1) ? VSZ-1 : first+csz-1; | 
|---|
| 135 | TVector<r_8> sv = DATA1(Range(first, last)); | 
|---|
| 136 | vth1.push_back(new MTRndTest(sv,true,NBUFF) ); | 
|---|
| 137 | } | 
|---|
| 138 | cout<<"...tmtrnd/starting threads"<<endl; | 
|---|
| 139 | for(int kt=0; kt<NTH; kt++)  vth1[kt]->start(); | 
|---|
| 140 | cout<<"...tmtrnd/waiting for all  threads to finish"<<endl; | 
|---|
| 141 | sleep(2); | 
|---|
| 142 | for(int kt=0; kt<NTH; kt++)  vth1[kt]->join(); | 
|---|
| 143 | cout<<"...tmtrnd/deleting thread objects"<<endl; | 
|---|
| 144 | for(int kt=0; kt<NTH; kt++)  delete vth1[kt]; | 
|---|
| 145 | PrtTim(" "); | 
|---|
| 146 | } | 
|---|
| 147 |  | 
|---|
| 148 | //--- | 
|---|
| 149 | //--- Remplissage par threads avec ThSDR48RandGen SANS thread safe | 
|---|
| 150 | //--- | 
|---|
| 151 | TVector<r_8> DATA2(VSZ); | 
|---|
| 152 | { | 
|---|
| 153 | cout<<"\n>>>>>>> Remplissage par threads avec ThSDR48RandGen SANS thread safe"<<endl; | 
|---|
| 154 | rgen_nothsf.SetSeed(seed_16v); | 
|---|
| 155 | vector<MTRndTest *> vth2; | 
|---|
| 156 | cout << "......tmtrnd/creating threads " << endl; | 
|---|
| 157 | for(int kt=0; kt<NTH; kt++) { | 
|---|
| 158 | first = kt*csz; | 
|---|
| 159 | last = (kt == NTH-1) ? VSZ-1 : first+csz-1; | 
|---|
| 160 | TVector<r_8> sv = DATA2(Range(first, last)); | 
|---|
| 161 | vth2.push_back(new MTRndTest(sv,false) ); | 
|---|
| 162 | } | 
|---|
| 163 | cout<<"......tmtrnd/starting threads"<<endl; | 
|---|
| 164 | for(int kt=0; kt<NTH; kt++)  vth2[kt]->start(); | 
|---|
| 165 | cout<<"......tmtrnd/waiting for all  threads to finish"<<endl; | 
|---|
| 166 | sleep(2); | 
|---|
| 167 | for(int kt=0; kt<NTH; kt++)  vth2[kt]->join(); | 
|---|
| 168 | cout<<"......tmtrnd/deleting thread objects "<<endl; | 
|---|
| 169 | for(int kt=0; kt<NTH; kt++)  delete vth2[kt]; | 
|---|
| 170 | PrtTim(" "); | 
|---|
| 171 | } | 
|---|
| 172 |  | 
|---|
| 173 | //--- | 
|---|
| 174 | //--- remplissage en mono-thread avec ThSDR48RandGen | 
|---|
| 175 | //--- | 
|---|
| 176 | TVector<r_8> DATA0R(VSZ); | 
|---|
| 177 | { | 
|---|
| 178 | cout<<"\n>>>>>>> Remplissage en mono-thread avec ThSDR48RandGen"<<endl; | 
|---|
| 179 | rgen_nothsf.SetSeed(seed_16v); | 
|---|
| 180 | for(int i=0;i<VSZ;i++) DATA0R(i) = rgen_nothsf.Flat01(); | 
|---|
| 181 | PrtTim(" "); | 
|---|
| 182 | } | 
|---|
| 183 |  | 
|---|
| 184 | //--- | 
|---|
| 185 | //--- checking | 
|---|
| 186 | //--- | 
|---|
| 187 | cout<<"\n>>>>>>> Sorting"<<endl; | 
|---|
| 188 | qsort(DATA0.Data(),(size_t)VSZ,sizeof(double),qSort_Dble); | 
|---|
| 189 | qsort(DATA1.Data(),(size_t)VSZ,sizeof(double),qSort_Dble); | 
|---|
| 190 | qsort(DATA2.Data(),(size_t)VSZ,sizeof(double),qSort_Dble); | 
|---|
| 191 | qsort(DATA0R.Data(),(size_t)VSZ,sizeof(double),qSort_Dble); | 
|---|
| 192 | PrtTim(" "); | 
|---|
| 193 |  | 
|---|
| 194 | cout<<"\n>>>>>>> Checking differences"<<endl; | 
|---|
| 195 | double vmax00=0., vmax01=0., vmax02=0., vmax12=0.; | 
|---|
| 196 | for(int i=0;i<VSZ;i++) { | 
|---|
| 197 | double v; | 
|---|
| 198 | v = fabs(DATA0R(i)-DATA0(i)); | 
|---|
| 199 | if(v>vmax00) vmax00 = v; | 
|---|
| 200 | v = fabs(DATA1(i)-DATA0(i)); | 
|---|
| 201 | if(v>vmax01) vmax01 = v; | 
|---|
| 202 | v = fabs(DATA2(i)-DATA0(i)); | 
|---|
| 203 | if(v>vmax02) vmax02 = v; | 
|---|
| 204 | v = fabs(DATA2(i)-DATA1(i)); | 
|---|
| 205 | if(v>vmax12) vmax12 = v; | 
|---|
| 206 | } | 
|---|
| 207 | cout<<"vmax00 = "<<vmax00<<endl | 
|---|
| 208 | <<"vmax01 = "<<vmax01<<endl | 
|---|
| 209 | <<"vmax02 = "<<vmax02<<endl | 
|---|
| 210 | <<"vmax12 = "<<vmax12<<endl; | 
|---|
| 211 |  | 
|---|
| 212 |  | 
|---|
| 213 | cout<<"\n>>>>>>> Checking redondances"<<endl; | 
|---|
| 214 | cout<<"DATA0:"<<endl; Check_Redondances(DATA0); | 
|---|
| 215 | cout<<"DATA1:"<<endl; Check_Redondances(DATA1); | 
|---|
| 216 | cout<<"DATA2:"<<endl; Check_Redondances(DATA2); | 
|---|
| 217 | cout<<"DATA0R:"<<endl; Check_Redondances(DATA0R); | 
|---|
| 218 |  | 
|---|
| 219 | //--- persisting | 
|---|
| 220 | cout<<"\n>>>>>>> Writing to ppf"<<endl; | 
|---|
| 221 | string tag; | 
|---|
| 222 | POutPersist pos("tmtrnd3.ppf"); | 
|---|
| 223 | tag = "d0r"; pos.PutObject(DATA0R,tag); | 
|---|
| 224 | tag = "d0"; pos.PutObject(DATA0,tag); | 
|---|
| 225 | tag = "d1"; pos.PutObject(DATA1,tag); | 
|---|
| 226 | tag = "d2"; pos.PutObject(DATA2,tag); | 
|---|
| 227 |  | 
|---|
| 228 | double eps = 1.e-25; | 
|---|
| 229 | if(vmax00>eps || vmax01>eps) rc = 1; | 
|---|
| 230 |  | 
|---|
| 231 | } | 
|---|
| 232 | catch (PThrowable exc) { | 
|---|
| 233 | cerr<<"zthr: catched Exception "<<exc.Msg()<<endl; | 
|---|
| 234 | rc = 77; | 
|---|
| 235 | } | 
|---|
| 236 | catch (...) { | 
|---|
| 237 | cerr<<"catched unknown (...) exception (lpk.cc)"<<endl; | 
|---|
| 238 | rc = 78; | 
|---|
| 239 | } | 
|---|
| 240 |  | 
|---|
| 241 | cout<<"\n>>>>>>>  tmtrnd3/END"<<endl; | 
|---|
| 242 |  | 
|---|
| 243 | return(rc); | 
|---|
| 244 | } | 
|---|
| 245 |  | 
|---|
| 246 | //--------------------------------------------------------------------- | 
|---|
| 247 | long Check_Redondances(TVector<r_8>& v) | 
|---|
| 248 | // Recherche des aleatoire tires 2 fois | 
|---|
| 249 | // Le vecteur "v" doit deja etre trie par valeurs CROISSANTES | 
|---|
| 250 | { | 
|---|
| 251 | int k=0; | 
|---|
| 252 | double eps=1.e-25; | 
|---|
| 253 | TVector<r_8> dd(v,false); | 
|---|
| 254 | for(int i=dd.Size()-1;i>=1;i--) dd(i) -= dd(i-1); | 
|---|
| 255 | dd(0) = 0.; | 
|---|
| 256 | qsort(dd.Data(),(size_t)dd.Size(),sizeof(double),qSort_Dble); | 
|---|
| 257 | for(int i=1;i<dd.Size();i++) if(fabs(dd(i))<eps) k++; | 
|---|
| 258 | cout<<"Number of redondances pour v<"<<eps<<" : k="<<k<<endl; | 
|---|
| 259 | return k; | 
|---|
| 260 | } | 
|---|
| 261 |  | 
|---|
| 262 | /* | 
|---|
| 263 | delobjs * | 
|---|
| 264 | openppf tmtrnd3.ppf | 
|---|
| 265 |  | 
|---|
| 266 | #------ les tableaux sont ranges par valeurs croissantes | 
|---|
| 267 | disp d0 | 
|---|
| 268 | disp d0r "same" | 
|---|
| 269 | disp d1 "same red" | 
|---|
| 270 | disp d2 "same blue" | 
|---|
| 271 |  | 
|---|
| 272 |  | 
|---|
| 273 | #------ | 
|---|
| 274 | newh1d hd0 0. 1. 1000 | 
|---|
| 275 | n/proj hd0 d0.val | 
|---|
| 276 | newh1d hd1 0. 1. 1000 | 
|---|
| 277 | n/proj hd1 d1.val | 
|---|
| 278 | newh1d hd2 0. 1. 1000 | 
|---|
| 279 | n/proj hd2 d2.val | 
|---|
| 280 |  | 
|---|
| 281 | disp hd0 | 
|---|
| 282 | disp hd1 "same red" | 
|---|
| 283 | disp hd2 "same blue" | 
|---|
| 284 |  | 
|---|
| 285 | #------ tableaux des differences / a la reference | 
|---|
| 286 | del d00 d10 d20 d21 | 
|---|
| 287 | c++exec \ | 
|---|
| 288 | TVector<r_8> d00(d0r,false); d00 -= d0; \ | 
|---|
| 289 | TVector<r_8> d10(d1,false); d10 -= d0; \ | 
|---|
| 290 | TVector<r_8> d20(d2,false); d20 -= d0; \ | 
|---|
| 291 | TVector<r_8> d21(d2,false); d21 -= d1; \ | 
|---|
| 292 | KeepObj(d00); KeepObj(d10); KeepObj(d20); KeepObj(d21); \ | 
|---|
| 293 | cout<<"-- Fin"<<endl; | 
|---|
| 294 |  | 
|---|
| 295 | # Doivent avoir min=max=0 | 
|---|
| 296 | echo d00 min =  ${d00.min} , max =  ${d00.max} | 
|---|
| 297 | echo d10 min =  ${d10.min} , max =  ${d10.max} | 
|---|
| 298 | # Doivent avoir min!=0 et max!=0 | 
|---|
| 299 | echo d20 min =  ${d20.min} , max =  ${d20.max} | 
|---|
| 300 | echo d21 min =  ${d21.min} , max =  ${d21.max} | 
|---|
| 301 |  | 
|---|
| 302 | disp d00 | 
|---|
| 303 | disp d10 | 
|---|
| 304 |  | 
|---|
| 305 | disp d20 | 
|---|
| 306 | disp d21 "same red" | 
|---|
| 307 |  | 
|---|
| 308 | #------ tableaux des redondances | 
|---|
| 309 | del dd1 dd2 | 
|---|
| 310 | c++exec \ | 
|---|
| 311 | TVector<r_8> dd1(d1,false); \ | 
|---|
| 312 | for(int i=dd1.Size()-1;i>=1;i--) dd1(i) -= dd1(i-1); dd1(0) = 0.; \ | 
|---|
| 313 | TVector<r_8> dd2(d2,false); \ | 
|---|
| 314 | for(int i=dd2.Size()-1;i>=1;i--) dd2(i) -= dd2(i-1); dd2(0) = 0.; \ | 
|---|
| 315 | qsort(dd1.Data(),(size_t)dd1.Size(),sizeof(double),qSort_Dble); \ | 
|---|
| 316 | qsort(dd2.Data(),(size_t)dd2.Size(),sizeof(double),qSort_Dble); \ | 
|---|
| 317 | double eps=1.e-25; \ | 
|---|
| 318 | int k1=0; for(int i=1;i<dd1.Size();i++) if(fabs(dd1(i))<eps) k1++; \ | 
|---|
| 319 | int k2=0; for(int i=1;i<dd2.Size();i++) if(fabs(dd2(i))<eps) k2++; \ | 
|---|
| 320 | cout<<"Number of v<"<<eps<<" : k1="<<k1<<" , k2="<<k2<<endl; \ | 
|---|
| 321 | KeepObj(dd1); KeepObj(dd2); \ | 
|---|
| 322 | cout<<"-- Fin"<<endl; | 
|---|
| 323 |  | 
|---|
| 324 | # Regarder pour les petites separations | 
|---|
| 325 | # dd1 ne doit pas avoir de redondance | 
|---|
| 326 | disp dd1 | 
|---|
| 327 | # dd2 doit avoir des redondances | 
|---|
| 328 | disp dd2 | 
|---|
| 329 |  | 
|---|
| 330 | */ | 
|---|