From e053a726eeee9e8f90b3290daf0ab4fa1b9eb3da Mon Sep 17 00:00:00 2001 From: NunoSempere Date: Sun, 23 Jul 2023 19:11:25 +0200 Subject: [PATCH] add example of getting confidence interval & misc changes --- README.md | 22 +++++--- examples/01_one_sample/example | Bin 22288 -> 22416 bytes examples/02_many_samples/example | Bin 22328 -> 22464 bytes examples/03_gcc_nested_function/example | Bin 22352 -> 22480 bytes examples/04_sample_from_cdf_simple/example | Bin 22432 -> 22568 bytes examples/04_sample_from_cdf_simple/example.c | 4 +- examples/05_sample_from_cdf_beta/example | Bin 26552 -> 26648 bytes examples/06_gamma_beta/example | Bin 22192 -> 22328 bytes examples/06_gamma_beta/example.c | 7 --- examples/07_ci_beta/example | Bin 0 -> 22320 bytes examples/07_ci_beta/example.c | 21 ++++++++ examples/07_ci_beta/makefile | 53 +++++++++++++++++++ makefile | 1 + squiggle.c | 43 +++++++++++++++ squiggle.h | 7 +++ 15 files changed, 142 insertions(+), 16 deletions(-) create mode 100755 examples/07_ci_beta/example create mode 100644 examples/07_ci_beta/example.c create mode 100644 examples/07_ci_beta/makefile diff --git a/README.md b/README.md index 45a3e34..d241d17 100644 --- a/README.md +++ b/README.md @@ -11,7 +11,8 @@ A self-contained C99 library that provides a subset of [Squiggle](https://www.sq - Because it can fit in my head - Because if you can implement something in C, you can implement it anywhere else - Because it can be made faster if need be - - e.g., with a multi-threading library like OpenMP, or by adding more algorithmic complexity + - e.g., with a multi-threading library like OpenMP, + - or by implementing faster but more complex algorithms - or more simply, by inlining the sampling functions (adding an `inline` directive before their function declaration) - **Because there are few abstractions between it and machine code** (C => assembly => machine code with gcc, or C => machine code, with tcc), leading to fewer errors beyond the programmer's control. @@ -184,16 +185,11 @@ int main(){ ## To do list -- [ ] Test summary statistics for each of the distributions. +- [ ] Pontificate about lognormal tests - [ ] Have some more complicated & realistic example - [ ] Add summarization functions: 90% ci (or all c.i.?) - [ ] Systematize references - [ ] Publish online -- [ ] Add efficient sampling from a beta distribution - - https://dl.acm.org/doi/10.1145/358407.358414 - - https://link.springer.com/article/10.1007/bf02293108 - - https://stats.stackexchange.com/questions/502146/how-does-numpy-generate-samples-from-a-beta-distribution - - https://github.com/numpy/numpy/blob/5cae51e794d69dd553104099305e9f92db237c53/numpy/random/src/distributions/distributions.c - [ ] Support all distribution functions in - [ ] Support all distribution functions in , and do so efficiently @@ -224,3 +220,15 @@ int main(){ - https://dl.acm.org/doi/pdf/10.1145/358407.358414 - [x] Explain correlated samples - [-] ~~Add tests in Stan?~~ +- [x] Test summary statistics for each of the distributions. + - [x] For uniform + - [x] For normal + - [x] For lognormal + - [x] For lognormal (to syntax) + - [x] For beta distribution +- [x] Clarify gamma/standard gamma +- [x] Add efficient sampling from a beta distribution + - https://dl.acm.org/doi/10.1145/358407.358414 + - https://link.springer.com/article/10.1007/bf02293108 + - https://stats.stackexchange.com/questions/502146/how-does-numpy-generate-samples-from-a-beta-distribution + - https://github.com/numpy/numpy/blob/5cae51e794d69dd553104099305e9f92db237c53/numpy/random/src/distributions/distributions.c diff --git a/examples/01_one_sample/example b/examples/01_one_sample/example index a8bdd5dfb6cefd749275808afe0677dba0506192..a22d80e5e7668542c56c350d77c79e5a771a5d0b 100755 GIT binary patch delta 5291 zcmaJ_4_H*!m49!B!5RL{n}Hc-2E-W@bP&x9Mx!JmGl=l0R0n|+8^x%Bz}7@C80>!8 z1&jY?6@xvS4cRnZv%B45o2DAH#nh&wCN9b5%NjTKBU!CzP4kHVQ4=)k+uwbEENj2r z`+e{I&i$Qp&pYqjd(V9{-EZ>UzvR2!_MsF`X(1K3dmDCr8c29AL+!km zK`H!0;lX&dt0QMq`#gZDLWOES^zHwR#|>(KkT&uzdqlH}XjT#0!+&9MFehjZdd09f z>;Li;(>%qBKc*LqerGBi=|OSaM#XU|9W+ki z|CP=fUHoV4PPRBV-=|1^NIMkiIHc8z%>ca~7vO6Z`8DPK*UEj3a^Fh%rT|~ZWNMCc zn-1a`6KS8pLAy+gjrALv>5R!kcTHucmgm;4Yo+<|bLsEm6AUv7XnVY;oilJoj4*O> zoC#OF(i?}nfn#?L_c$Q0T%3~U5ZLwUO_-i>+!@F$_gi=U1M(ivUC{bbr60$Ag+6O+ zd^1p+C(hB-l+59Yl{~V;NY$b7zQa}zdw_vU)Mi7UNi*|?7Gs%36ScY1`3`Csx80*D z6^q8j3Lsd7L3hWwXqq`GdBYesj;5Pkv<&RFF>E|72fG7o&luKBU0@G^?H$7=(QmIw5y6 zC~kdTjo4YZiiu+M^VTGBdtn1QQqx^omY%{RkC~R!i>U<#)5a_Q)C>r= zrY+?2sUqz)p6JW8NeRMjb>LbL8%RM{(zu8R_IbBytJ_DbT@`#6?R7QMV4B@f<)#a1 z6BGV(R2>s=-ydVcBuE!V*_b4;^IiY&4cM@F*SbQ6Z;#SGS1LI?sl1c?9VugJDx3WaU2s5UWX z6NB!EH2Z5csF&%UC)3dmy^;YN%fNGRAT6Rmx-S{Q`JNd54vaA$w8k|dei2K0Hw`({ z6PG}ge({wo54X{Qbg$dyMACoIfr>gJEK5Fz9(D)~G?-@e}x;kS%k`gdgBEQgbjO)E<@=VhO?=DZxFno;os1;s;?ORnZ5T zdG@ET<2|-0Hb|RU`e|SMgxC3;OC!S6%bULvc3wn^I`8qqrb`f_Bgpu9F*^G3m)Mc- z(udx>_JU)=rd`m9J$PK)zmEuY@4>$H`^?sV71t@v-tm=q_*|UG-`19RgJ3g@oiSU3 zul3uX9uOxFDeUOBI5-eqs0SiL)33cc8Q(GA5>WUU%0 zIDb{i7(^@sreI!CwW!;-({S$!`Ld#hf3ecfvlbcpF4H$zPZ*9}23bfiO>9d2LlVaw zyn&+h4!eobGtv&SPx_W~BY5eAsq>cQ1nwb6f zhBp(bKl_Eilh8e%f|s9F2ML?~y0EnsB`n+@&%Adfzm&jn`6nu(U!VlOnn0^_mhwLO zZO(**gE)RLE_qJ^eV+4c!+?c;Jh?w{9$3}iTv|EBHGL7hVE=J%1J$q*ael0#E4daY zQbm_(K0tXB>6S!4nR4K{*qtL_)7>`66qy>9lw-r z_-o_3zmR3x?XMcQ2gk1T$t}9=U(!DGCs}4EL-_G0D zGmq7Qsu32xyCw7^R%o%Gl32HWXMay3>$7c;qr`$;W4_L|8mm^eCzdw|e>)If7G5q! zk3fut>#(AIqUCslC_Ece7m4x$->K*-yq?8z;m6CPzYo`xM{kC!!_hl6V%PNZpsOl! z@3N#xFP_KZJABvJ>dpUeq3_YDTjJgdPUG9zCmhVU(LX>9Svn)j&A?IMHsIipEMEoo zpOxhxu=kuS?`D1Ae&7Ib04M=Z0nLAuWjD&x9N;uyBQON)0agRe2&0kdKbPfA2Ckj! zh0p^+LKgM|{gN!RZ_fc_@&@Z;^=9MiY6GwYC;@AM{>!qw2G|SiVEh$X{t?hUf*Ak~ z0<+`LM+vM2j-tR00DEx^Tm;tMl4YAofC0{%2u{Ooz#ht-F@^7-${8+G?`ga_Ju{=I zehB%|8~Xs;!8*tX7P!(Fw;R_T_VYA4_EwI`;4GXP7_BM*=Z=Y zy1R_UR^N`ea;txfsl*!Wj4!c<3?0^B*y;~keWg}+v0$~8Tg@dHGZzaoc7xI%DK>3_ zDD|y;q%Gc}f4E~k7JB(+YAN*E`}5c{I8I+bA$lD;kq#FI_+GjO*+uqw0eeLIj79>K z7NN>{Ub~;=2#tNcBOOLZ`Lt}_>E)d`t#^qxYQGx zF;=0XJ#c4(t7onD77lMlxzv!`z{mBAzwjK;vJt4}U@51ISpC9n(zjbua zZJJGtcEKh@d*+w%pU}1W0i-Ug#O2X{*Xs5s6@kT2S@P#(>e_p>tHf*An@4Yzl=4^U zR>>5;m9iED%C}=VrXZB#+J;pY7DoM8MV992P4G6%r7fbnR`CFBT=3_))^+dE z;S$%%vzjmVQOugM?SL0F{y07@SyQ$u@RrVNY16Bg0ZW6$%CQoJ#;P=HPHHK2d2%#u z5Z{cfDOChGhiUn8#fH^!YrLQu8F zI#=4yG%q&Ii@x=Lt@Bx`9Xm}}JErlz2Xg#-jqlcStk=x53eQnnsold46c%wv+eLcm ziKC68%kIWDW?J^sfwehwu*}<@tFawZRaV*L5Lp>fd`(jMJ{7nnx)&`DeWz~Fo;{@P zhJa>Yr`bRGK-^6#pTb47#OU>}l`2J-hn;!s@wh2%P;J((2}j~d+whV zZJ#!FNLv#6p>>F9XMM!R8vb%oex{9jQd?m9;V}$enP`&{bDF}Yyxz$qM7Mi_8kU}D zuf|6-et{Nsp2mkJtHdCZ#oA*Ue^rZAPi852_Sa3SX%(|zl4YG{8|(YDd@Ha0hci^+ zptgm{!s@lpG_PD&zv8L2>sPI*f3B}{;S0R`>H5|c1;G`o*FM|O^i=(`tLs-ZJ===I{K-lGGLK3 zR-!$Xvn_hbY^34J+4qaFUX&x0zoZ@IWZ#ki`KoMmbxF|fv#VZp9HRWHK#D%+Iv()Q zx~hi~^s@nV^=DOEcn9TI2kd%v)ve;lrI!Lw}n?2Mu+Dk_!y#E7kJ#PsB delta 4802 zcmZ`-3vg7`8NPS38+P;B%_f`O&1*O0K?0FYAR)Yx?81iJK!8AqE{_qgB#;6UGzhj} zf(W>yk%Zp|&R`or7Oe;($c=)Kh1x))x8J#sG{xzex%WH& zcmDI=|D69n|Gm5Iuk-D{=i43eAuA_0uL4|6N8ayDUXx08eCCPTnXCvqlgY}D_ccW- zqK2GBk$C_c5=s?WnESujt5f7Y+RoeK1FBX))e6vI{)$P@jG*c1Mcw?2|7WR3wbVnV zh~Q9{o*|V}>_rUQ zYj!gGIxdnuG#{51JtO+8btWiSqsB|)S5U{Or<`gI$tq3r70poqp5_|uRBeojZ5_et zX&czBVB1Hqk@OnaBVap6utquu_6x8_N3byzYqC>}DJC{Ff{mk6u+vR3raq0OkIJWZ zx0+fbDWtRQwTw1SE$6tfls%2cT0SpA=UN%ZnO%O4yJ=3p1G(-ldG1|u1v%o6klvbN zYLg_X!JKO-rd(@*E-{^!Sj$}raUAz9+ru7q6WWR&N#i&%7OFz)P2=i>9i@-6t`NN393wncx(kB%=#P>luEHb#Y}iJx zCroixr^@z}T`*f>{gA_Ei=l`z5|bRU`KrjUf0Ti3Su1SP^2B|*>{Pm)n8iz!Zri}` zp=WGaQFs0=Nn!J1!!J{)JlUD3-~0${`Lm3So1KE^+lHv$*1&J0nxr}WcG{ox68|+# zwP!}nzoX22*U4s5t34*N(}@XdX|H`f|KVQ8n(|li;0E*c(_z6p_I=p3z2acLE=E5- zCJS3+e)52#5@ISL4813sZ!30N9L#QYjCuEN=wP-!i#fFnCd4t&K|p(9Jap5%#Ge&| zq2svxS&pzd;Loyz%}YXJ>mY);<76&WLM#My%$^=u$&#U+3X+neXQL}_gHrzns&J0x zD`<_=rOO_qpflV07z=k{A6!+oBetGG&hqkucduzEX?RW5Sr&Owhkv%H)&Ei!^>->u zPD|XQEJ${9v7s9V#K&kuvfCHpAxyWDzE5?e?7GEfVjHllEH}hWGGR%}PN0mmbW~?= z+QRfVP+6dbJwODGi*t~-P?tqZ+H-U?%^zumS@B`IlQu4X_f4#0m#l*+M)6elr1S%P z_NB0S>^D1qFh8>&TeRaR-u%Q%=mf(E^go2);F)g`^%FGY8b^Cv=9b*vwplC>PK3?( ze;~B^PZ`awC!o(9vDk$jr;x_W&cL1bu$YiA<`yg2=PNs`NvxDY+@9Noyw ziFzzXUW)C(SaN4==TFefSRckotG>>Z=8#@O$v zf)%sN%e@U$+=j=6q3X8Svv{Law;4_VihiS#o>`z0M%2iruTb9*?9W<_) z%6P-i`2gdu!DX?BYZ_{?{?oU}7aXb++V05*U0a>^BYUR!@#XFf<2v+sp?LW>3iA0D_Pyf5 zN(R4_q{--?JugWsfFf`Uu;+p#T>$R;N|M~bCof7;JClL00uKS%o1z!^0Z;@QarPWv zOHwLO02TmyfPSEO3HE`3%aYWr<65|V-@qUnz*3F?_w`FsA8_!hB(V<*;Tj$c>;V=5 z2Y~``3({d3aLsiX0`5a5hJb;vBn>echp|8pJ(QSEU@q#}N#I$a2oz8&3k)=MLJr?V zD<{|uJ?F4e+BTuQzKfJIJ(0a{1sa*!>4&E$9Hf;MRf526952TX#At&R7{O z+gJv!6;W6_Lca;z7vNZgp*|M>x8Hun6+Tc|TR!4Jg(1*YkAyf>xia5&- z-JH0Bi0OP(oTY+2k$q^fu!ktSXk5z>dN&Zi3b0Ijm##9-y*AF`i!;iuv4n7cD@pZi z9E~lRCC659>`Rh0#<4^>1nyCp;kZO#09@jCl5~l2y7%MUb#X=k-6rr?!RPBZ`fHId z@)+d4Ym{6($uNlSBC0F)81j)2D=1jJHgZ3DePK$Ulx%WiHp`f0lr||P$`9!hB=IIa zK52$AA1Pu%3ekW!LHYb@H7zh{-;gV)s?^1Aq{dQL{FbqpiKazepcgS#LvNLOP+0oW zx00#M6CY5&=-`yxCR0_JE8d7|!5aJ4WKEz)%jWaFbh^yNzfApQtN8#eC~u{pFPd6S zS!DEzusYziV>H8Oi!gD$K)X}W3@B#gdjuEF8EJv9lIN%kx@OJ-wDT93q>tyejHQdd zSBXm@qv7feA*!frJxZh4kD~!WY(+qaAk(D2-vn^^c z<+O*&XDZ){uL;)Va|FDi@<-L*7Fvb4q3LJJ>YP^%dIb4q*_{cHn2#4&UZcsU7?{IG z`4V<7yjRufWW^axPJ1A$ofejBH8MSXTTf?Tc`HWabxzKB*JzOVO8d%v+H=#DY1)xH zsPgS@g~w5n+nXxinXAa!UY=BWuj;>6tv;4lcTpYX)N`R7_AncbkHisUYY*Wwk~M+O zRJvO18arBHyA|NFRsM#$zuFmVP@kdw^!KRxg{r>gF3+Z`yj$JPM^*hw@D?tab_?;& zT2<#lyfS!J)d^6SV0S)*(Rg^?TcFVDvE??&qr23(A=OXqt+0oU79U`95f#f_QF%?B zWl^(AtJufDyQxc;z<(Ljwyure2lprM!hD zCatIm)VQ!@xKwH-b24=;Y(X(Ato4`>8>LZWZ3(TYweWq^RqL^8!PM|M2Zd|LM`<-1 z6>#jLUA&&U7J1^eI;f$S1{Qf}MO^{a)LBegJ*%O{x)PICt?;K$UCD5158$UcYlmo{ cPU3ODx5Ok?DO?p5-dAEOoTo_1*aO#p0HHn1LI3~& diff --git a/examples/02_many_samples/example b/examples/02_many_samples/example index 1e63c346debbd2f4694f8c4f38d8ce0cf2ddb3ec..627577c8221e00ca7c68a305f82eaf4888a1bf62 100755 GIT binary patch delta 5382 zcmZ`-3s_V~mcF;4+irxu-E>2D^TvQYbfiI4n22;EkBhbsUItc;Is!tLMKlnTth?Z# zXeZ)8&g@1foBf>GmzhszNH7kWO-!pKA+zhOnn}jY=xU6KjWH21MoiSZf8B?}xVzQg zy{FE9PMxYcr|MMQ9{L3z_$A-xjO?>AoEl3)*&)XHU%i7Jl<^TEIuAQ9Wl)a>lE&? zitX(Gd&JAYjj4R6Uj{=`Xk7%2OhSCzCTY+{>LJjmc zKo0`#n}iza5YTf#_fJB@={(SIQ+UL{B-Bh1CI`(28kmGyr~qgS&{LC8+u$bCO2Z45 zX~xV#jvJR!XHcf)s=djK!HM!eTztOZ2$6st9VB~?kUX|BVzdbr#fU02#=HR%jV_Qtw0nicQ#lVTi##ut* zcGeOWJoL7GFP}=aQK>p9hB~8K`5)0&Q7K`W-$~NAP^W()hJG32vW0GW3I>Z)j0{^8 zgYlM-7`=f%L<7-_`E7JJ`eogxPI}2PEv)m7((-L5Ye^R!;Zyd*qrsQ%&^%WX*_=Mz zy$M?9+(5e=R^2~2>f@Zt%nL5aH7+rd+Rm~*l5^g z<9E{H`{RT*hS2f1SkJuhRbjl%W!$3-FYH<(6LIP@?P)UCc0 z3qmOD+NbCQm`(t)3AX0wpBhDS>w7rU=HAvgiA zRl^0sqpyxFC|INX2D9-g6FSMa!RUc1w@!$)B;+T@ zr4Ld_cR!4pQnzLegmK(ck&w>$1!~ypH-(+;NMPO}c4Q%`eTMTX`wUW^*gNM76WONM z>H4lx>-6|ADddJ7WjbZ?cVK#2X4& z9F~=<8->3e@vic&7X8Oj4R{~NB2N-cry50}HK42;c><65dJ7&+<~Z+)a{ou(nsWaw zZ?)Hdw?^#!-b}z%zKna!|G^8NA{xi;g|sJT|4|scn)y-4YwnqRC%ZyIjhpyyh*6)o zAW3bY#*31)i}9dipiP$~$qiZyGy542dKmO+&=JtRpyxn`Ks`uTLC|VY=VeK10u@2G zf)0W9F*zoEP{(z0weWNRz+Qyu9O&7gBu#)meN~d;a5+1R73~8}|4NctK=*=nF&y;A zptaW}=>WsAxrRYUkasMh9JdQ8%m+Gxq}T;|0LQ@?Xf2K#J63}6TS>~(3s6AiX#jaR z!2*<#lg7KLBFAAka9)x)o?3GjM6r9HR@ilJH=Y4HmSbCR@DoYGO;0XNav0A6U@irc zFMlHppM?3qpJVb16ghioWgA4vj74E4Jn-UfU>P!^VAa5w*~wTFuwh_7Q1x(Ta$A8V zU7~-S9qrnSYA%BCmy^T?f%VbJ*^Uwxvcm{P9w;Z8-MR~5W@GtOE(z>qU}q5?W`Qd) zTe@|N&F+xde3;qdH5+|XX96}NLe-R!J42s`+IFhTja^-f+F?Yk0diK@WLkh(zL2C` z8Bs8Hn#~9HjEb>={}|W-V2`psL*jXJQk|J?EfLQEWdBJgbIbI1Q5zuFoOy;kL^}66 zt(Y^*z!rt)20b_D3H>=#2Pl4SY*-T7;d0N-|C$o!#w^0!ic4XO;cH3ycP1(|TlVOR z%+B7BB6HI2&~kHnkG|OK?lu&gJ-RNl+iOnunv+V+&LXp=+-xjHlX|olZ&A_QBI7|+ z+pwhq2DY>;Lp+Mf(^V#6@@!_%i!t4L>?AzcpMNU zd!uHXIHT4~F(SW1xG-+l|A|lKV9HxY1K#KmfpgQVUPnw4j`D{~2xsM*X3`~Zsmqh| zu#kOfMR9rLEpm{jC_2OeqnnBx7^T1HN3MYS^{lmGLC5zGH%8eD9oTW@3$t9>*N&#M zMb&Aee(2=UD+|l`yL4k=7Mz$|>~Lvc#G3wxa)phdqBtu@bE_t`U~%ZUaAgf@AgC2f zE)pYIJCj~3F6H;ro#HgUhmuRmvR_t3U73phfZA+8#V_KD#Y(b5xBzV{aq#cav0?`e zlo%hZ09DNWo=Mi-R9kbn7PBJnGthflZ@9>^BJWMmsaP0nfL&^DIDh40SMh-qh4-in zD2=?O4%a+YdJ%U}}>9WMmM^)6Vwv&?qq$@OgkI2g0phCEM&5G)K zZ7)Bkj!yeiC`;9EQ}rzmjry+@ycFqK(bSS(4+Ivh1j|Ftd_gsoq;A7zHGD5JoK28_ zm0lGc`okI@|Bl+CPaUduEd5r~;4hN3+~K;W>RHrfrX3L9k+Y! z8c5B>*(%Pdc!@go0u>*dt^h2stSnLS{SUdYmQJI;b|z)>F>P1%+6Hfw@8i`!d%TP7 zahunyuWwvazp>*7Ee+iIjepg)ZgYdmJ~?=3@mAisp`m?EzI)C3jjfH%^$o4-8`d

Rj{x2FVUAmbce0|w^o}VCBMV2+7R|cxlII07Bj{QT2qR~v+ zSCJK_~PR+x6AAPJwT>f!BFSDj_mQn5w@w68kLr2VCYbzeP? wMylPES7RX~p4uAsgEL1gRC*qweKnHK7)C!?G1nw6Qoz0R@rrqD)jXE=KMvuMiU0rr delta 4942 zcmZ`-4OEm>8ou8!12e*)TgQAWCiAMReifFo`VNz>iVJVvw z_HKq7b#m8j$D4?d#3y^cIy=VI6W-b!hNchPqp&VOTt-`o;yL)(^lQI7MO@mHHbOV~>1TNNWPz(R^pxq8NKvb{t8&1JZ;Y)DvJ6{==VS)!S*^ z3_X?ThYJg2{wbMnGuVYWvVKt3yJ5FcZx{Y5bAOPZ%jD-qxqmL~XLHH08#ATQPF)o3 zF-(hJx_lA4Z^g=nTKUPK*0TELtJ%WF`dZ2eoI+~@Z*LmQO3TZN$Ir>l%+0?gj=Cj{ z%pAp6uVVpt1~5G{AkqdP*MsoX$xlTVp!03>w0Q7a6}uSO3}Po$u4VHUq{+a}tNMn` zTbe-V*)DpnwF10yPH`A;%VC3qPoQJLd{#6q@kT{uz9>cTs;Sm^ ztI8XF^j6R#4)q|)wsq5#Awd9vmKkld)EE-lIDpmDb6`8bZXCb{(h;zS!L|=zjr0rH z8(=#JupuqYz^3h16V7yg56{aq22&vm*3PihJrenNYTZ2Fm|O4!xlsEwi@Ph7~9_v z?S0&k$mK`<4JPT7#kr2BWuJ<$A1mAE<?pN#kN%<|RK%QlDj(;jLKuA=VM8UwuD(mAH(Y`*SSj zf0Wu|76})gvP~9_Qn~Fp-AD(Wx22fhy{-(jIQT%ywuc1v+F}0%72B8T3a-;m`=Y|_ zn7P*S!=XN}CH(_9j=qc-!ZqaZ5Eo)(*|<$H@pBVDy7r1q{fg(0ZQOHth-JeD*kH8Y z#aHz*T!_!W2AMwy=l;u<*(I(JbnU^f#AWZZRF}A{eU{mNvG+&hw|T!CCVrj+T53-W zEaj)BjZVeH22X-0UVzcL<#g4NBD_KNSf}tbjfou^@lgyik<*QpD&YuygX&`oJ%1a# zrE=I|V!owvmsQ@=OTXFEzpp8KN{hjs;-j8aHY6`|mLUKh#EayNONx+J?S(N`oy1pt z7wwEm7ADfCF=@fYe6;UY3Ez}U7ChVk+_!B}iJ(v(@dNV@qG(1DkC)i`J zUFHK31ca2I4DpOiSQ8$i?gXbUHIgnROiwh+V|WA*{ypO4U-%%atXdOrn-WU`t6?I3 zLpu^P!hh_;CLWV*FvKY4cYT|4#RwLE5vA`{$0%pX5@GUzVToHh+(^B(LwoUUfT0WG%e+iy8FWDipY!yn$3lVZ&U zDt48{Ohy)?`k~rl8ZF8lXjs?qmiGGJZc0Hnt455lt(ya|9F?CH=Hu~c;!2 zywwm1kas5@qR-s_ydw{W=Mr#e%B>uuEnpxgj(44|p(Vb;2$BL@TTCNsV7PZ`{^0tm> z09{c&>=(X#(RrrJIlL!;zlwy{vWE*z{L9ISt@M8-Nn;?3Cnaevupih6Jbp@&P654N zOHvYW%r}zM#&zI!U^Q?Tun~9|*bdCbWjqKh0rmr{fpw=PX(_P%j3jLY8nN;Y9cyA< z>}wbA5GCm|p!YjT>IZh>2)OZ`(hu|k>n=)C4X_hf4|M+^N$Y`joXQTM7q>|S8gVgQ zdg!?7)j&6j*$dreALK|_q>i@J2;0z$=^8#A#WeRl|%m$zuvDXW#rUlbI~TW`rZNd?cf&BrIEG?uS431 z_2ZPrZ=%MQ3hFo#(uaQDMO%XgKgjLI;$keEh# zXjI%Rp4(2OY9ac0H(|3Y@|UnpfZqcW3Z*Fs9|uLz~g(rJV&?LMQb= zwvjL{E8M4k*kNsQ^HI?_jK+QEjem0TR!1wwO-mMWNqNKlUd~Y)gTpH8e(D*wNPvgp zNb>9SG`Fz6>v*AEz=W4Pwj94YoDbDn05dVot;TT-U}j93TofG;%-l4+s8q1hPG}aE z1hmC6ym1D%02a>N8T570L`Qb+;7;x`nt8A`8=_ED935a|K3W0IM=ixqI$jv~#f`!I zaO!XXx+!P84U?CR&vIzr?3(Y@s*R5}!^T4|jV~4cM3=^AVTR-qn?w6^PP6})EbvuK zFUiW)BCqw{L^yQ}uMB?b2DES-;s}s5BcU>abV9 zTakLJZ-5_&Z_v&Oza9;(?>*X6Vw=~a`tqwMM!q{RIj{2B_&DQ@*HhlEtNdJk_|foc zk(&Hr=~-WIEgUcxsa ztfwBYP1wL?N4siSuTH7GJGRM|T~t(NbG)nC`0)kIn^wd>q7`K}LFAV6eo}3c)Fs{u z8*<14)}5$~ouB}-s(k%Wg*c9o@Rp?V>js4^gM>*oQ>ogU)HQa}R$&PLSA)kh;6Dua z&0AK#bm7W5zgti{XTkDS^Or8{3MyYCbUi;MQJ~PNL0z9tT_SXOd~*fi6m6TH72au3 zrq$3!Jz%x!rlFgXXJnbR+9n3l>=~m?b!tJaqpdUYZ&g{X+WM$xMw8G>+bW#2uEJ`v zhAR_UX=_D3wN+TLfS!u12rWYzraCC4GSjS;eAK|%m5&K}>Z!~M*D9fgD5_Oi)K-;E z>#D3Kt)h8pYgN8Ut5-hyyej|Jo55G8tm0KlzEjfqO_VZow24hp$h&CX%zV>bQxz>1 Hk#PPOkJ-E!pQH0zLCp`Gyop);N$ zgbX*ErkZNh9FPT&5d8BB0&?Dx|4mQp0`fW9CfF>JY9*;w61^;3G+Mb6G%LNRTbA_y zc{-tbIzh&uwbu7*3P+zl^{cvz&z+>s+K9ThPYs`!<_+~UYNGSu9_%z>y{}*S#1I(7;NM*goU?FL4#=t(5*nbr=TJ9 zS3nN{?VExc=oO%6fbN@uhS6D|VaBlV!6~ST!i_ds2-G(PHB$l57N8?j(8!^7-<)jzKaF8>dA=^rBI2p`f&gf;dyw(DV@ge;%`?Ypup zJ=>CO@bv!ZdOY2YNR8deq5QL*=cV3)<1qE#yd%pV??wMvaz&=-zQ03_k&fy&Vi=Qp z!;2JB(1bp6)Xas9-RMKN3NTWyzL=Bz$&M)LS^eFI!0+&N80B-Kqnl^$m`Y)9SL|~m z=#P={!YP`JoRM}4_N;SGwbWDa8fPW%R~=!}vjuNJ@PBYymi-%W%lGuJP?0r1e|wzb zFAy=&9KIFuglN`tI##ad@$80^d-@ht@1cHcbgBr`kYTu!8^p{bTZ~yg&=cYdzGDpy z8KPr5?+FhsI&2*fQm7^>MJLBnS5&L;*Ys^va%lRWWO+iY)jt+Xzl*g;25)%+4oi{^ z99tZV@s^Sl-6*_5gV7HQ+vs-mFLhtW&`Y+NpVN4LEwCq2 zWQ<35XObR`X{6mYi|#izrBl&lpi`Gk>lCuihEA(#zC9^S>N%b|wgD!5z>XM~?w#-9 zAd-&TBZVDw+U^kQY0^F`VsA8pnl&;ZS^^^t?<8kzk?WV3YJ+(7g9%AY{{yBG6$U#@ z1iZzQY(@&j-7f}geB8zd-BJHPe-{|i#V8CZ-r5BNB?M+3f@4t_+#wo`O$^825?7DC z4{IcZ;$vonEag$}qrT|a=_|qd6QCR0LiKSDTN8G>_;fdF`Gi;(o((VXsn4dD;${o! zbS!SB@O%0?&S821cK($&10#5jBI8rSIsoC;n3=zsO5$fD>MikSg-TkS5NF&Kp(yo@ z)S3_z`{{Mw5kbIMdD8svDu6j*3GGX8=x#>Pk%Sd9KDY2bJOUG9%@O~@*b*>TM03JB z!QQ_g+SPdhQD1_A6k# z6NvW($vb)S5;pB|+L}10EB`%ldp~s20B#rj{)G6*yRh$girc1+<9|lAw|*xb`65^n zeoz;9gJ?ELJw9`TSoa?dV(ka${#k<SRyyx8Ap09Ab4=HK{G#-2c??%nX}1x0h5&qA>}K2`?3 z)N?roierdc6+%!ou0#wh)&N&vZ$KrIDlp zDt(~Yni0m>n>SFf4)KfR9ra(nOh+H6h(0+T5np&kmX&JeUkT2C_A(_Urx;g3h3Q2V z(z4{Ibu^tmPA1*mP-;rqnmHKC*b@=Z=DrFJTm8DYs{;khJY{Ogis4K=oo&>NPjcy_d*OUnf>O{y?9au{Oth;e-q-aU-Kf`TbJC1 zFc_Ma!Y{TXdkZu$70e>b=jjdKk28{|S5F|y?EJM`Vj_dyt zU&HY$xTk;kB9*1B(Y+i(|CsjUV0KZKX?a@Zv|Ffw{!ebvZ_`$GE&p1U`J=_X%3bZQ zs_@=E%yPq6wdB27iP4(?)c{i6q$3dIEkYy)m4czSGGU)4|PlBER9RNK8It=PUx$=WnfyP{r z`i`5{~J^%wq(;3jyep#Lbee#kl$K#N88Y|iZn)+;*j&d!PoVCYgBjb666OItfuh(Ax*v~&anKq(Zmd`dh96`(PcOm% z%##4}@C5TwTJ~(gNh`8#A^XqDG83pZyD*BM|Fp)Bb$jq1q$Ampg$K^dGERI-g24tv;(sJhL^T2JV+MKvGHQ-)H*6N|>jZdZp znE4x7&QXv6#-C;jfjy;SeBjRm+Yjs^e$Su;!IW5Q;#*6?e-OHVrw?<=^|!$dl0A1p zNFFksbCp)*W`yuX;kri8F+X-O&Q>g$FNPbk{I6V9 zW-{;96`NvugNjXwdx9%WsXO&0CTCAbiOHqwHaXpGFO-kCFrsN{iT~! zJipj*0Bjr96h4s{!|LX6tM-ANS%4xx>j3mldv9w;9DGe%Lf1+hx^;6YqO?rdN42H1g|}#1X?fP?=oW$7bf*VQ2h|=| zRQ#MewImh4L+49vLNJ{y!|QXI;ofS{Ze-<}3MZ1zy)qn;c~kZr=z}WWjYBPO%02~^ zRD3X5ftgEPG^@Ec`_EMEIh+%DOHroE*e|#i2X$`qs^&RSF=5J%13OGYxlR8rw25?3 z7b$!~S>?8_^N=kVA#(<-vGA2Pu6pq~0FjU!97N&k6wAho}}hSz!})bDzcm)doAQA|)#~mu-IZ@wS%wwOJYSvUBq^>G4Nf zp8_yHM*|w_w{SUYsC;p^kgI&mak=h&$H~y#Y+Ad1L&Mq)8#{m6QqR_Je7x<^P4%j{ zdFax^n}wLh`i`~v&b8||wl*|xsBc|gzqYxxqkhv9LqAy}2toLYL>`ZQsB!5AVQ6I8 zIzgDEtQDCSpFS`njmA?u&@=oi9GXVc>EMdYP_1mD%wAeC->6l~-ITU+!M%#CRpusY zU)d!v%Br%Fp)!WnRhr3J>9jcUp2V9*$Em$CGeS$eMoVMqwaSc8?ID4(dZ}`!V5Rn| zOp8`rHR_{-RhdTZV^JQBRV|=%RZd!7Z8mD8ou8Q12g>1@H;d74d5UnD1U-z1q=u>*r=mOV%hq0645ZRFm<;=!&0^= z^xcfzLwj`B?O?}K5w&UCsqtr1%T{zdt!?>Nnx7d8_!Ir;x6l24jOgq=XTIlspZnhT z-uu4yzTXGEuL-@c3hlPY#waE?w*stv#^#KdO>0bwx@Pld*}(ZYXaV#`-Hep~=gF^n@{>tF zKRWN`qM{F14aHBndSGVZm4pAtX}Nxo9@D#G_|CLNF}7ZI97TKd>B1{?R&N#l%b#Qk zuu=Q004fg{D?A|cTV%f7U=!xc`XO2Gh25qAoA6he`!N9|VMap3iY4rkl@B-7$xjBgjy0}W%@(g|tfS1pY4k+koz0V2rMIeV^8CE4yu#aJ z{%vV|Rz6?7j_L8#vjAp5qzyo>2jQucpNcF%=iBCKapSitZYi=E#7?VR>*lRUlYw1Q z^^KdiHiHOYH$d>l9jpYZd1cW{mxl4w3h>G~#i9PT!zMeQK*vJ(tY})|jf%>A(TdOquk>t`rKf&c<7H(IH|7#g-_7#lz@f$aplX&4(wAA>ytwqqD; zq?=%W0^2o=4W)4=D~&aUhV354hSU9E=Yc&mjI~f3*k?_l)UP-8_)YFW!efgi)tFbp z*r1d)o-!=o)l6qZweWOQ5t==XWBbN%9xTax5?ron0SwvCHQq@*ZU+*VFaI=Sks zr>*E*0%M{9gtzsA*hdSZ(sV($Xk(P4W~-GkPg|HnhN60yNHJx1F?OvT!6t0 zm^9QB<9WuA#O0^_O(yA#*|CnNXTOTDpCsECbGBoTlPfXYjHgC27z+*YHjpU5T?t1+0D)zRI0O z&iyG4^FKu$u}g%$E!Mk*W8}5Iq#JLiOV-rj_rF&Ln(cfb<=8?4`)#oQoXTv=bw%G( zw{1zu4$NF@UtaPGfqzPH>Nimtu*l{;;N=Gt?VWFog0D^z7_ahsCL z18ZR-enqb&Wk!5Ah)p~#+hB-MEbRF@X}geiZqS_Z^;0*@&vha`Ew=>ohL<3;4kF8+ zds>H%eS?(mr(MaJRFZ6Nb{%S)!VkgzLGy2pc-qPj8qJOku$RMN^|TzsA)a>(%h`vW z`aP2luDSrr19j%C6a455nr9wf9S>>me(s+;x^$XLxqexG(bMuR&XT;okNVZ^7?Bs7 zGA5)hOPjBq{3*`~r)Wm%?{$xa)90zv5y^~E^K>s?CSPh1{d3gP=vPA-+j9-2=so`M zw4N5{^--C#Dt0D1#9e)oB-bf%A0%;FA8mK0nLZ7X_Y(6S9d$mbiwdRcG}1j0LSouu z**WNY6z9cKkOisEoPE{2u?eNAbT@zbpj<5uM%6oxV&-Wpcrch~^jsagM3+Y=25-hY z2m`~`22)J>+d>C@oW5szI7r1`C?xnBG;jP1FOO~WoQ(*DvAzm^%?>wv2dHcXyU6u> z+rnIUPkGx6Q2=>&;$iy4_4m69V0bYRhgPogKRddE^M+sZwVdz6oT5vl560ZD8xurn z8FK;}FGv#gj;Rh>iAd4Yf0iSj?y{awxX6)#!Bz0YW=jdrmFRN>FmiH{4cq3=5(9FM}9N5axDM^|LSv)OC3xGqwHNX>R zBIk=36faSm;U@frzoFp{>JI+heCZG{3 z@6@qo=E1&p^A1sxJ^^~Zk)$DD7mk1n-z!5v53v5SB+UbM0ULp?E0VMUXv3-O1bT3r zM4%BD(-{CASG^YKLXkWKJOTU`xCtfOf;Yktu+~6#Sy#sUUAGB z1Vn563+7q;w$ZV?sN&AEl7x>Nx$&Mi@#_F*S19>==Mdx_!#FSWFYxO>MXBSbFJFMJ zF~iK=h}vT#xY<0C;1;QNc&D&>aM!@40huQ}S$8-%+*lbRPq_h{5s`jyn0*Jh2WVit zwc<@kn-G7GGBQW^Ek=4nSTd}8!M%%|bn(85a7&Bs?(iJ_T|#iUrGmR33tGObG%Y{V z;6!59(Bt{>b9rvNkgCP#=Zy!0=Q-e#B<0JD!ttEvfLp6_D7lOm>rQazIH&tCJh>*^ z=z-V)zKO2qR~pVis-@C`!oXyt^oz^%L_t>IY=~k%y-~2tZ~$U0rA%-bc*@T~iVrS_ zhfpNA@fy7{q0BfNV#hVynE`i{_64@kqZHWl`kU1<-<$#-0-;HUSe3xUIQh)dutI0_b@=0d?iN~y5gYxyW0kKHeObrdnhT9t$>bDw_o zB*}LQhAx%Q!ABc!6Xhi{HbLW^3h!1=Nfq~I$D2y_){aja-sF1^Y~a1Znq0%a=otJjr`??x`t&U9Oy%Aouw~-e?nopQGo4ISD^iGB3tDjRPFsER`96G_o_S8YSVh~78dub+_x$AxJT+$TB;^+3tz0T ziz4uRO_u%qS#TT!cz{P9qEqU`ethflrj_eH@VKnHC?SZpS6Qpn(ouQ(WI$qx%m^>BMX-|HY}c>n>9YKAXn3tFMMb% zCki#9Zt-fa=l486b)8`0|9)`0fBvrnrFf%zD!dO1bbQ*ao-@-M1=KM^kI%P25YABd z%~7BEnq9jRqY!!$eP&CUwe?k?)zve{1y0W>f>J3_018rD;8 zO*U<-$)WZdi%F|d9@<$`XwoW>k1o~}-hM0iN|aT+PPz9;I=_kX=1efLsS3G*R?aCj Mt(dN8afpQDe>D>RSO5S3 diff --git a/examples/04_sample_from_cdf_simple/example b/examples/04_sample_from_cdf_simple/example index a11109c117625a008b347da1c2f015c75734e4ab..5dbcfc6765cfb49cc0a0062b80e78e2607a31fd6 100755 GIT binary patch delta 5926 zcmZ`-4OkS{xt^J2ahIQ+Sr*t`epduVP*_w{H0Um%F1kh(Ni>Sas6<+$7(q>*G>{O9 zdlduu($u#46l>a$r|A_#dW&g(7Jm|RoAeT^KBm`3Of(MuR7Ekys&n5nGi#{#-t#;& z-#PF1bI$p`@BGZP?d97$_!dV@jg3=_MM< zlhzDTPSV;uNH`vMq1JZK{#V*7XlB6_LUm{BBjEg}q-jI{4*k^C*20mdh_t>u;*>=hb(O`rb;* zjJf<0R;RWwhw&gfOrw25EOi^_7*=koAu&9e^1^2u*R5asbS?cP+#NN6o8d1lDqcD! zXG{U@56_6@1kS*PLDRq)@r2_&gsz3#>RF3T$DmG$x!i;2<<{5X@*!Lw>I@60Z6kL> zZ`UPMBa&g9f;KBmzZO&v85dh!9nGQ^2|qDRlQ*dHtX=D5`T`f1rw1k!R6Pf{vRFu= zOw_{|!lhVxXUM-fSqMU9!5pnN1cawTlbu$ZqM}y~LJjm3(49ax4MM|d7tqguZW)A{ z=yRZVfbJNCM$lcLqs&p!?SoJYjW*k9G0?6-s6>@Ow*fse2(?kO+0Msy?lTVye^p8| zrA_6yekEfVWl6VVSBV&Vf8VuxiqdjjOt}rl2eiQYSv@MhM?#CKB^;+Se)B+k*8}aJ z4sLrZe9iebq_wZPa3hwRcYZ~=ukn&eenrLd9U5yjSn_ikG);r58B|}5#z)Ne@*G7` zVCx9t1zV%Q^?kz0$Q@tcurF{xKF!nu7a($26LB&FpWa8WXsC^$ooW~P+md0N#Xm17 zb1*Oe67&Cvi+PqS4IGxgLZjLjIOuO~W008Yf&DAK#)Brg7Kz+1T8nS_3%;qlAg_dE zm3Sqn+QI-DVAnrV?)QmV!7gUMUtD550#NM$M&#>vF`33gro7;yO7CKC<-sb-wk1|H zC14KC(V|KU+hJCTD|(3IZnR*^)u{NIjZZLAaH}@Xx7C=!>YIbLW~En5Z9qrdL7hUr zK$XuQL%*~o2?LX9zb!dyU@~OtvMYRzg-2P{7dTfN<=a|#0uA~5J8;1Y^!c&z^U1Hq zj!$k%VneU&KB(pdX8389n0sWSrC*eye?;Z61&by^Lf#6!esO_um@2-|m<`G%cd+Vy z+IZ)e2JRs4E~pwK!C&Fq8kUi%j-Vvnkr(s@&dbGCj@ufRzB!1Qolm6mu_G&ICF&9C zLc}=l;y%`w>Sll7Yk2xKcnZO*3!&;*^eIA$P`#XzJGrNAd5p71V! zP|NM!XwC~VPm4bbA~A~Mf5YdK9G}69WU;Rme&wWk`>=>#e9yeT(D=drsDchMAapgNr_rL5+muj3*LptnGiM&Oh zB)Djg)6YLeSDmZK;j{{!@v6?=o0`rm@%MFV;2rtTH=3^G_JwAbkoM&G1s zP8-cnv69;_d2w9%JwkR(p=>_UWDj33)DV`MD+&RF@Wb2W0zzHmvVrF|4H6dI+d~keDC-icMV3p|we_XpGObsBZ$K zu|nu=XG-|P%-b#GawJ5~Lsf3Rp(uSDX&@=pz5{m$@x=z@T)((DdJ?R_rg0KwCy(GC zqw?f5;S(EeOir=92|@Y5b!`a0raj3S(Jui)-o6aG zTiP`o)Hum0Y{Fbk^G22osY%fM4EIbOy zzoDPMSBW{_{0Ke#wO^cjbYlT(2e6}IiL+iW{+SJ>Zv%pLzr!nMRn8$c3sve+-*oV? z)i}U72;a_!3l9B5bs4Hi1IA!haiURA+IHcq%k=oL0B@iR!{!LD9I=ksP=G26pie*~=F=E&FM1#>}`g5%Boi7xnE(K7g7$n>sRvo14JF?Bb@=>5Lco zNP1%UYr+OI-5UN}ZWk^c8)CrB{}0r##;=Pn*J7LZc7-$J&Cvrm>++8J1Lv``{L@5x zMl9qxx;rAp>;(gQ(LAMRb_p#*=&Q`GIj`f&p_zLPHJIxcU;Z`B`I;A}*RW%D)x$8$ zPOkSrwP|1m*p!;1yKukqHyihWs56uIe@o{^{$cVYNcJVGK8A;Wur!t7#{Z31GW;?+ zl{4vh))HaAk@B-24XZh)DAb!(KE!Nnn;i>$i99h-7WHEY5wy#)4&J!W#m}3~5nimH7M-b|B+MfojmZ z1R3t~nR}{y;+mkAH0lJ{nbMf5ohc{06{Udz?|g6I4sJW6YXO)0N8V-rmvg^HG`jB^ zYE$$6TkL#x)NtV)_gTK4eGa5y4xX z2WKupga5M_oowHHCRCnT^I0i`!u2Et$m@m)GJbvPeZR zrxiAg=b$Gc4nqfFEOw)SjU0r%0nGQWSOUi`3@ZRWp$jsPdej2+p&qpL#W?-6b;dF$>R2s>-&>RQeCIoi}<4v)1=V9zF+SJ$K5D%b82@AAC}c?) zB~M5+yP#Eu)+yF1WuapC?NvCs|AaaV)pnuC;%GJ$SzOz~N-dd9EL4r*B^Hm+U~waM zh~acKVx<;S2~^&NzD(W;#il!`NytgM=CO4K%}%2=^r0uQ1VZ6kciWJU=#-Wfrdl02 ztj=+`QK`jEt%a%lMA}!F%b%n(s1GH6O0LzTe`>)MYH1+vlvL{~Jz?0VAuH{)Zps`! zle(v*^0(>ClqY!)d5dfm6vF87shfCL=ND5`c^JOtu{&?GpAHWc183$kX{guEx8ogW zh~!q$T<;9(@WvZN&O`gW_5@GN14?KbI4kGKp|8Euo!r>}iyg_jJM2)7&5yZ`1oK&T}q&$ zj=S|-P?sO@WF2q8eU+7djd1Ic+sw}6!!^L8=f))3J>Bk{8PdYVo0SZ6F|yVEosb>X49`l2*s5_k#0F zExk-09L}w^v2@gOMVGJ99UIyn?u77Rn%pR~S&7vAki|7*Yd|PaiMoTX>VXRFn!~Ak zhTU19bK6F1TA`h(lyQev>f9Q9ez6h?#+r~#8tEHR-=w#C^%PMz5J>M1+0IavxlWC6 zzYk48&-q*(->b_99~kL-I^LnXBD5n+pEh4^#oFnU#NYk<1Tw=+RUHv z`tQMYvtpB;Ub=im)zTGf>z-J>l3Tv^shVX^uhg6CJ5S9z!Y4ervUX{Kd+GAEYpSYO ztX#8v<*(?P@Qbl?_rvvkXL8wcUXTQOac++Fq){KYPKQ$W z++1sDOAXOc#Lvq$XX@K?CKb$^U~bo2?exsN3HNvK&`#b!-Sg^siM-{x)KM-`L%BP` zgO4v(JXV@iKAd{XrI=8#L-cqe3G;IzLZ=C~;z{$H_$2C{pBr<3L)JHD$w)#)jyd#4 z(oF>w6UbfRrnU;n9NIfQw7+7){Zq#?RU1K&_(v3BVgw0~jyG?bsewny`{)GoO#OgL IfFn}>2lH!MZU6uP delta 5380 zcmZ`-3v^WFwf@hUOqj`I9!VxMlh;f_LPA2543Gz20y88eCoq_p5Z*{qF+ic$=p{w# zBhiR3Ks0n?rPp>f+S&zPY}3ZB75cD)ARtR_qsv~gUaOY zv-kh){qKMOd+)Q)IkRs-=>L__=SeC}X6guPgpD-pF8qacz^Y9LcE4K^_gaQFc`bvI zgyPWNN!rw)%%RN%kXZaGtF(FLqJO2QL~TAoj|%Q2N!OBeEol+$Zc3%7Fl!76U)bh?a>>!qp7Q zfSH)6UNxTkZBl*b0S@99+5JNVIF_9U&#`Ul+{{Mxd1u-N^Np_1VyU<6G6E^}Ry~u&SWENogY)c6?m%5aa$S@3LL#r$peZV;+QFd( z;Qz$qDz_?%0{#7({vrsaV@yUtNRij1$g6UcuZnyOm0xHoxAFBaqKtHSSc7MBc#pb8 z*-M$Q4wwa^vJ>{@12+|=K60ep?t5I#mlPS5D^Y|!mxRvC@{^FZ+ef475^Qu3%fnMg z#wB}UbcmZ7x8G?Nxcwv0M(6gPH$=V!kym4Q4c3o$DyXU9-GbB zl|D@OgOPg8&vP>C+Gdxu(|jlAw?;Qxl@YtIoBMQFr_g^u)i1k`PA6xGWwYsea%N!} zFxFYrBK1@yaaM|axH&=ES(S>3JoPJ_aV_>dYEGO@!IbjMqZ!=m)rSuE@ucE#ueaMv zj&HGz+Z~B75~a*qc|SDd-+!qnJk(~OX?#xLoM8wD_A>3FS$wrLFme-|DR{1CfA~^ zLZyTwABB;kw_I?kCfBhw+PQ+K`9pbt5>pGMO!$Q(R>I70VR@yRVB|5hI&^1f&Bddz zrw!;$fD*gxq)$}98)OuHt#dAgT_cP8YbuI$^N(^R5r z5kvHUbUDNku< zQaxR6o&gK;BlPDCf0F9>@hht1e(v~6YV%|ZwbbLukDtdKJ_+H)bu>TI=a%tSV1KL| zl{;=OMJ(hY9H+~fd5LcWQX^)Oy;mqN zYlwypq;8*2_#^fEntUsq zdeP{0^0C1m`kxc>&X3y* zKiYcL{^XAlqMnvX~7`^Weski2s z#ps>b+AZ}QLou&<6aKu4vqq&^$2u=S6CJkr6fdB0`yH=uNe6cbtz`5xH+1PN=K@ic ze+C@|6{Y&Oa#TN$1vTZ-*2K3L8~gQ9$K)r)|NQI+xe?)Kv?I4cG<`~^b5{!Gq)eVe zop~FQ{}#{K!7J$0ukec$IV+!`&+`_z9tAGX{zUIPauYcDJ!+gX#X1EN9K`hzbxnCp z_#YDU_lT)+v^W2ul4Ws>eGld0C|d#!ZS%5y+h#QM(2&55_a;7#cZ9%O;m8HFkY}w_ zU9eo}qn{LHSOxJ5cQyI03p(Z{jk%R)l-a)gbE0yM4d%x3AM*i<)LMeyAw?nGsTZuV3rE{}*bh zk9-kY9EyC|BK3aj4B+PQ)UUYv@>3rS{`p%0ad6DqD*UGOywJtpU!lN)e^Zno_~2PZ zSr6I>+Rf)9I6ZEy#&e2Nz~}EP%3jb`&;y|Tpi$6a&{LotA1KO~oSs*d`?0PcbQkC_ zXg_EW9={CQf#VE|tc#6+IS;@uD+*o#S;s|1$pH0VfXb=VXbB==! zfwD10c>{C|^gO8lvZ9pYZD9m-6R01Ztsj)3%{vgf0nnMCtv3{9nHgCg2Geh*8`Gx< ze#)NVHeWjfPpE3fT-O+u8H2xa?!vF1cF#yIK7g>4=_3A$`7(Z|fE98~vN_sp7Ul>0 z27Y5B^wx|VVH15Zqq=JvfZ2%8*n5f+56Wt7j$LAd&GUpQWb^mLh@&=7waqaP>(#-x z-Nf|~|EU_Jz72N+OF+ygVtarU0Q6&k&b!FM+wZJUYTb923ad4GL_T3Qj!B4yw*dAarRYb!M0E+^9 zK*xC4qrj%2QJm+Pc*5pyu~{Ur{PQh_%;j~lx!6}?4XCN8vz3Zl4DZz=7C>h+1uD|5HBbs;#^_*0wPg>u0zA;@^I$S=s2`c* z9NZ@ie0$thl_hl2@+zNmh=%|JuM+C0q)3%d2-5zl5+O?Og3qAqRV7YI{|12T>R3;+ z=J=eFUI+Y}g^%sDVNQc!p@VaL!WnvZ&N@M&(A=j6PtElRF!y=T?LEUkY;PGRW<|KO zLvCR^t}atNljx369V-6UnAuq=s5{M%$K9>N#GGts8oeK?^^Q%yRmgD1x>zNJ=DBgA z&^(vP4a4_hCQ)SGQ{GnngWXuis|jChtd9cKZV_QxP+jaazT6C@EflJDQ-8HZsHFYX zb>f~18m}&f{p=dI*Z2fCl>eeK+=n}A{ux6<8lc@ZK5<|g{idc?V067^uJ8=anqLz5 zZ`%~sju?mJ(3Ug)&^b3 zrNSzH4Lqpp1oetEdeu%HAIVmis$b`vdzyqgxA}QUQ4wOY@M{XxxxKF#JLy(@+9LM9 z)K=$id(V*Swr1+mk9B-FM+2r|&c}5fA5ohmIbtA5&Uywvhxw6qj**1Ns=x-_WY zMn;>c;L_}cx|HGTfDg4rfv(c=0SZa^ZJP}(J)VQAE@RK?9vS~DmFW0mx_;v=9`)*Y ztA1jm744;SlH1FlN<3?R{Vr15J1p*3&v6)1C*<1$+9vva&5m-Z@D9HQE^;05o7MM) zZn;sfAER-{&~sRuMS%tGMf{n^cOUR;_Ki5*rsJ(TK3@~FS{-Kv8o*m7A4`De@EXXu zPuJ+yv#+O(J*wjase|8H@Qg4S{~WakgLUEhdDUyC7f-KX|JC;0j!k#px8|<3?Q7QF z)A{X9cMq0_j|=3z{pMhP!z04rzQ$HTye3dwQ?WB@*7h}M2Guu}IE@}=&_aqdl~^r$ z7qrl!rV8tJeYTy1=8CTyrO`0gQ+;z6z0>HR2b%44q1kEm@d#-4{bUT8imfODV}M4ys(@ zAj^_ctI@E~rkB7=y~#`XTZfM!iY!sYF)Qs@T5c^}pn*@*>q{%*Ub|hJr3rMRG3y_A CNw+Tm diff --git a/examples/04_sample_from_cdf_simple/example.c b/examples/04_sample_from_cdf_simple/example.c index 7734adb..fc8e219 100644 --- a/examples/04_sample_from_cdf_simple/example.c +++ b/examples/04_sample_from_cdf_simple/example.c @@ -87,9 +87,9 @@ int main() printf("\nGetting some samples from sample_unit_normal\n"); clock_t begin_2 = clock(); - + double* normal_samples = malloc(NUM_SAMPLES * sizeof(double)); for (int i = 0; i < NUM_SAMPLES; i++) { - double normal_sample = sample_unit_normal(seed); + normal_samples[i] = sample_unit_normal(seed); // printf("%f\n", normal_sample); } diff --git a/examples/05_sample_from_cdf_beta/example b/examples/05_sample_from_cdf_beta/example index 8a07292b7ad35d9147d03e8a680d254b9f2a0ce7..2381ed049e3ad2059ddadfdf47ac1984be22205e 100755 GIT binary patch delta 3939 zcmZWs3s{s@8b0TT83u&;2ap|Z!wdp4-ndw-VY$x8?I+PiMMcEK1k`jD6)){$VPWR< zh!6GcmiBbPUE9Ys*R*L{D!0aK>9$2`rKN#p>tJ3GwX`+=eb4{@__TGNXU_MX_q&|$ zeBXD@5A7fGjUBwzEmhV$J@~7LxB*NXaRVrZXN5Od^ksv;l$Y@EWUeohX_NGEf92W3 zn}=^3e|%|W+V<+tTW0)Mp3y<8LTdRgIvrx?=fq65(0CruI1jij8dqbAr*)w@{8Me! z5vuTQg5!^AYd7>Z&@9tTevtN#i9RSY ziAgL3QH8Kj{7jh5+B%F`KqYNvn1xgf4o;lHhF>{;%M#~>3;6{DYkvA!Gh<=whOvH9 zTf--egpE&Ihcc^?ihK^oD-2~2)_>l-D26~13_etEVN2ukcqWJ}HAY(QlF#ahY}`rAIt}T~&MKX{PfTw8)hZ8@ggSTok5T z1a>%svETu+B+uc)s3CbQ=d?TdMLwIV9BC1kuj+=Q(}f|mIU>Vc>CjKQM&GBn=%}NN zKSHjQxuiIvrL;7f?idj9TCZ-S4!C2}g`Nt2*(+=!mFE5K?z7O5ogV6Rc%)Iil$sJt zYn-vt^((Z=nS!opcMgm(!$`^P=#``ONTWtkk8||s7m?w5`SPJ&MIL+@$B0{>$`}p; zKDsR3raTdl*KN{u0zxMM+1{Z2iXPM)TIh20Z-QLQfK6ne01kq`p_g5*C4ZKC zXR|2H0b9LX{ZYLP>-6L44=@WJ0T&!cU4sTjy!k7wQ{vvl+c@|!KAoN! z{0He^9Jzl#r-v>QH!Z!XuH%EOV(lS{mcM$Z=v+bF#E9a@PTW`-WDNgK_lMWKv_$aJC?Ds z7E>&Uc8=hoGZdHg_j^Y|v)iHN6WLH&s?!^gIm_J90yq5v@5(tKud_ie=M=_Kp3-THJu^ogx-o7z<`7pb2c#XVlA&e&DV4&QN@N z!3Mk!pqX7tct^b=FAV7Yp&ejqS<5II!kBkTsqZuItZ-v))YqHWZYz^Q^WrKbc zNehnt3VpZd$B;&6?w7KETh6bouX_rmCh<)bAKk8ht7|;^Wu_8V~9^sj7E_)|^5>;-iV|G>Y{ZRow!*6LdeQ>uXg#E^yEbpk1K& zxH8u_s+tU1j-H^QJ{#hxJ!Zgou)Nn?u8K@Z%b zQN<;uGSLm>D;{B52dV0ST3NJlBE0Z!>SRh z8)%q!z_VhjZG$w%8s8E!#_C!hT55H#H5FQOn$3k)kF?sF*N;yiSgMI zdLnf~9HmUyE2L*86!RMYpoyKFFQD*>m-#sQYlTx{!|8BEG51jPC_j~LP}I$de!--zhB zO3~BoT%ygGW5t}8ffPqezwEtS;KJL^WZg~QBjo|^qFqufuM+m93ky|9{zdEBu z^5@n3nd9o8n3WxF3DX0x_@f^TlcWPqO03O_=(}Qg2PW38<;irSHan)KKv(UfHj_WG g?l$MB{hSfuE30(rUYa&%M8wFM`YMrARZYtO00_WwKmY&$ delta 3482 zcmZWs4Ny~87QXk9gb+er_%no`BtZBT39(wWc59FTLD8au0*gvrL0xw%u3CSFU9q$X zjwpKSE_STN9Xr+SIz{7l&}s#LMvGd}Y8SUv@MrZ|b;T%K?Y4b;?t3pC$KIK|@7(X4 zd+s^sp8FoO_Rx-Vv?Y~$GJ8|{Ro$SCNMq0jVKg&HtCgo(ZvoAxL3%1r^{|yY7xA|b zOZq*ho2DL{w&2!a#o5o`@E8As+JL$A6Sx*&qPN%+jH)Q`sH4E8ilVzEz7-a$^65oM z+oST(k0rK2dajVz8eFqN?V%^2OTA82wW6jPCTJ$Zc1@ms49PDn7@sqHc;@g?FrYE7 z=12g>q#|lW4L()UJdkJ+Iq6C8$+YsIXu+~D;;#sQ%c2(qF^jtjSHJ1c8bSiecgnm; znyU$ii`bgA(p*JCPZCr+;^vZMV(xE}X=x(A9Z^T+MPoT~W;pZh#~dN=xwoR&xkdPX@^ zZ}+_2$v5thO*)xLC#sGLeLmT-Bat~S)bnfCq7ANA*;!C}bRh&j62)ph^%Oeqe8iW# z?Rr#o;m>W4^6{l^yV1v&cM2|a<*7fOiYA@McxxT#Ef!6_h-+J*E;c%N0*b;+G`h10 zwpbImNFN-qTIeA3T9d<)6R^6@9t>43N61Vt+QyCh;Q5^X9{+Prg?6Kq)b3v-HR_?9 z)LPhJi;LJOCnSYgrv42&5N5#(HB5<(r+Kg>HYIpGi*^qx@63ad z@m5o>l@R{TI*x5$5pvKb6-j*%7RIO0cLCxPXg+)rZ#AT#M!4`m4j=(;#3$)hk}+0G zD?EjWgcSM^UQW0|Z$NZnOwf06QbOtn@JeD-^z{F*g;)bDDQgTNS0apYhoL6XO3y+| z;#A87X-O#NkJaxICjQ9Q7^SH(ZrEnKThoCS!V56LJ~VXReN3T2GQkyE;X`kOeFsgt z>EqLHz5QQ)!#1o={ZEvC>pc{je0auRdzzk{{}vDY7!)QCg}x+yb?T|+Jk|ooeEh64 zo@Vz+EpL4b?WH1^JoP8BgDcKsG)J*3Iz3}fF1wAE$E*3wUe-81e#+;|V^Hck#{4r! zcl(&+?3Cocc`?3;YGLvTg`! z*o7Gq1vTkw=%?UIe}_vR0*ym!GaH8x@&>kwabyEZ)x6^ozjo|}!R+GfrXNZ?&BJygJ*_b9N4BRhyU%;XnQibIbvC-s<>pqi&yEc<=|3xq zE|dqp5XJe3%1toBGz3J#Y2d;uSBsIu?_JDrvIlX+R)zB zD~gSXHHg~~+Ys9kjh997IMe@26#F@{nzWR+eMbU^~64w<`%s?EtDvCvjMb|M9 z#5Tk>MB@!n+|PK#Q;0o?-ynA0L_dfFw=fYEA?>$Cu@-UQ4sMLta}NVkW8iqZ+=y9t zP3Noe8{CFMyBf-KlBo^W1xx zK}L%99!6UzM{7f!9#3Hx(_w7&+(biGz{^x;Fy zKX~Fsx&m!F1&H>H*X~47Ky!Fl6pptbUJY53QfV!eO|o)A7Oa}|Cz=UQC-qT{_wr;1 zM>P;XB@!%Cjxz1_DS5QWn>MwNQWXRhzfZ@)?qVy~p9S5;c{B?`OH$}|7*pbKvMUmX zCs}7^Vc}~<2C`4#I~_71MUCPK*jr+vXW(J6iGBmuO01zmma;^F2TWHGjP;TH;C6x#6@DOAqRD7>mN*jqywY4QoU`#Fgibe^^Av4SjLaq}tIA-F?BRzeQ+?Vv&V^Dzf1LXQ|a( z?PnF2URg>03TsMD1xpoe-E;e{R`?brFMq$UQ}}*%VR5`G19=~=m72^SDcVX!`$t9V zh0rpSx`Sz@%3gqMoL}!w@@SdN5kcM z=KqV(icIgT6+5Umc-A`3+gjO1@!pM@;|S6O%3d|z#d88V?gI;~nwzQf-v+#d`{p)K vJ&dVxgcaq;rhe>Dc~@2aM5z{5&KnivDVMcd;K005x+}BfSu}-%m6rblX$OB0FFcu8VJq6x9;ec!oruZx)W z%(`cv{q24B*=L_~_PKX@pB4IlBXqk{x}8id{up5oE;&%}UzUDLjO{=0w?b=LR*Xr@ zq7-3f;IJvi>Xg}Iya2*Q<&`nM{2O1Ue>KGTA=)Ll9FnFbX<8EfO1Nmb8dG6c=}E(q z+<#~3lxFD^S>hT~GiKfQ@^AKDY;?LB&GKMj_0IP`=%Y-!alVnqqSC*52*YA_`6m97i#^GTHmkMdwESwMz^qr zw!qdgI$+$8x&FZxcK_xLY}=Nttzot1Y?!rdZlhFFopH;y&0AY(t?BFZYm?PbQbGT1 z@=Rt17Kb0h6VHrzOe&Xv+Q4|tP{-3nQ*Qf3w0Y$B9{3TWZepKmd_?8r*)8zA2w=78k7odTJb+LaQwGw>= z0!opSax1W4bmVZA%TniYb665LMdee7@2KdQHx{aRfBUcD?FJh}L+UbCb5b$Zng(}n%|ZK- zXsbmT(&Bhrk7rk{|k5xVoMm~K* zQRH7xWjmkNR98VfEwgX3BRt4Aa~+$2GzFv$0s0 z1h?MBBqNUUi}YuAZj#h_B7gVZE__ZqJvqWjI_R0~Oy{vL z9*Au$_?dh&od!K~tDZu*!s6xAqmo$gI%bjfL^|JdE75j@CfTtH;-gQ-bfR1*3fZf2 z!sjukHq+va?9_J1)d1Le02aZ3oI~vyzGQgk`_}LjC}a8PC0CAVF%Nk+6{Tk+E=5!R z?q`ZJ+(yPsugjN)!}$FU?Cnu;Rq_kqVAJ>l&B>f%IG0N6GIND0+LP(Ec_EUIT#I>- zL4BEdNjH$BFg8Ckr|Ep=6rqBWvfeSg?WB*gGOSg&5xHlVog`+vGG<-nU2(wRE?+e{ zPX+ARM)GHSh3`^R_R^fkQ+OX9f>H7AWAY;GXK1UeJ$pM1WY?H_p(0PDjGSqX*`Hz> zdsH2CG0SazD{`I|@-B^v1)n_fnYjOBM5yzYAnv||Mr0HLKQBeb20sjo`6npLH?6(m zHF5U=$fSc9F8=qU;+m7tZ{5jt^GERfL(@`LP><r{A=xHpTTkQKM8D)+ zRw|qPJ426!KAGH~m<}{%FO{B~;+nA-R&e^*3*0HyZFpQ94Rs|KA<>7rjF}+nOyr|7 zU7PaTIbVb3aJCv!1VrjuHK_%`GdEXp-1t>|HODW(W!X&^^4AE(R5SJN_@;}BLL>Rh z5_aKjL4NZKS~zuSd+!H|!Y{MHioo(fLojmV6;_tSmP?Uq%P@pp9)LU;c}=~i4s#J! zTDzPA%9~g^ek*;n{o@zlGj4#-U+&o>@EF}ajM48xXJA-&>;R;xHv~ zwtav9!9>pULvY~&_l4hHsG)#;VE@+Yu=u@GfmMOkQsh-Mqk%Qp&pydApk6<*Euv1~ zo1w1CpXM?aSP_i;Ij}q!xfW;$L~bmXx@PPHyeu^JCJ&libsk&q*v+_BZ}Hn=-wz6A z#UCkZ5ZbT8e_U#8;9W(@NBjoPDoP#b5NI>!hZLm))PGJ2lZc9lrhj=P`?rV zQ1n+d;UTInNzciDPf-dvUX!4{D^5bQ6i?4PFfGuIlDTOkT&xw`t85=0iB6YPwKt=4 zt|sBVhw$tMHgz1vgWPimW=61H0(PIKhpV4u0dphJd7w zJLB|20^V{S1;~ZyP-HlZ4#7mmjq%lc2CH$W%m>+S$oR?}1wxYPZ1gp&+ zu$e=N(*YlbEJTBMO*8%#t%pe{&0IZ()-xX}%6iCo;S)~6qM1KZ6do~NIOYRZ3T&r} zU{au^zS3pU?M8$Ue~o;t`*(Du>=M(hc8n2wdsHu0~f>6htN=}c1+*17x& zm6Uzk_%ND>>FF|Wdk?B>>X^7%nkhKf>!{;FLg8N#ymZsfxl0W81bTgLjxa>Sb2kVQ)mQBm z?0th(UI7-b_+1d>2GZ`bFuUF@)>*CNn3Z+^OFWU)QGH+m^#;=8M9ks2fGfj?>-)|U z#~jQ&oz4YnJyOY?LjFya#>y#Bjql3p^f(tz%9d&uyo^-;&=b`@yY*gtYZ$O;MyPn6 zOF-D>6?*g!EM2Ef(`lgt&?%>1&08S+CtaCWh#=$(_hq72^~~1zp%JR6Zr~a`s$uTvrof2w2YFgnLAFq z=DUQ~=~#`6dgq&O?*M(YzjI0OGp%=zwlUPx$W^cbxFbHii}F$qJzi`YUdIy@ieK9( zxwNI$Z9ijPT?msSnF#4 zC3pusVMVGk&$4QHUbE8Vjq&=)Kd$3i72)(xRuQ%u6cHsIAuF z;a35VMfbFrt=pN0o|9aTNt(2Y?;@0PI=aByUZ$Y~+G?q%AJRJA)YuCa zK3TUgguGiRYeRYer3GNN+TN^8i8rN>#6yCUQ{ie0+8#iuhUB7ilUo<2LBK|Y-`xk_Q z)$?-qF07wr z(KB5)Me1kY&Zc@+ZKlEcb^~ttrAud0bg7-XmlipS95K5(9ihRcg(>>#>GXn!QkKoM z>e(InaN)8&!c`hvhE%R)MV&U0qoL5E|D#k+!G_u7Zz!VP2D?SiUU;8qn0@;aknraW gFa6|Lt_X+8v0}DG4#u#pRKH@j;hx3xdfnvz2V#(%V*mgE delta 4936 zcmZ`-4OCRu5x#F>VVA$%1(sbFVHaExG~rJXtpdA%vM7cH#Hh6at)ezYR>f#id%$WG zVlbF7KQ=W7Pm>(e(=-s06k8kDpHyR-bI@vvNsJgZS&To3m|$Xj{qEcMjKS0E@y*OP zbLY;Txp&@My8bG(zbUkuBRx?}8V)~U{SR#)^SZWC>lYihpH2@sW$_EAEEFlEIUPZM zQK3xlmjw{6a^%YTnLnr9D!<%A8w69NTefn`R_>X!yF7-xLbvvU-;Q0NovO0L|L4vt zyYrGZpgOvGYK<^@!kxEP*`BCfthv_xpN|T*(fWWCBj2798DkqH*HN@LV6^ZioeMAt zXQbyVWC%>BIe}(jktDb9C$*~0LbW7+B+1=@W@>`X>k{@`>G_!SY?Aafu%E{{^**eV zL>;PVI-o8avAk{xTfK5seXaCl5NnyIZY^89&QnWcgWR+sXtHY(D=sZ5D6Gy(%gVXS zj=#%|OB)ZLDi(l$j0G|^a;*mCQZSw>=?R{KgA*IR1j@_Oy!{1~M9r?s@)k)BWVaym zAzX_I>LOYldJgy-zH!s%^$d$50$pG zah=ZE$TQ*Es?JQ{Cf8gW4{2&$x)oku7qU8_GZYS$PCCEli_+$6?>h=~%(np*QUs`r z>>A(IFm9*LuYtcza}UlNvtz3=He2q%Oit(@bF`@RMDvv6v4)=e;x z0(v9*&vGA~1D(;I#k4gSeHe)#+%GDp0PzJ$a)4tDXWmI<6 z1QQQFe18llxl3U!8VjZ{HrNVr85mcKdLbv1hU$$7CxP=%dF!=eubjt2N(px{+`a8mJ7ScZC;xHDx#{@UHfpr zr5ke$j*Y)H^U}@*ZwRwvuezq)eiJV@@j`XbcXr6{`C&8nTpFg^vIRDn9cbrKjfV^0 zcd)6|r5uHG?+x9&B6}EA2k^hhZXVFhDYEMabo0DEH_Tic52wJy%X3L#=7gYPB&D{M z8e=S>)4}>?!lrR;mv;5WeY?u}&cCa-XqzQI`dR;uq&VlQx59z%A=+=TMM`1s`Bn4(4mB%Cs z(`a?f=+FWl?k*_%mQZ$_)#NtfGXK62_hLXdKl~CTWJ!IA9*rBViZs%WxJ03j-ixzF zltJY?anm2cbUGh5D(nHt8u`-r`Y3wDXki&m8F5Ax8clKWmXJEUXs{)j|7~S#(eeEx4LbAM1nTn5p~(oJ8k{Gye@5Do2p{%|#)_3fi`Bc~fr zqr>>ZuA+tTX6H)+=X;YJ_T-Rl{@!Dp7ET`f2H|BIug;CR}^=Vz9j6IIiqR)O0 zT&*tshHgteI*hYR;NDxpub?+(o+xep8twB2WhKv3)qF$SlH)^4puoIvA(Y-t{+r4c zLMujhg}#6nq2J#Q%1kjmSOPCPtAqWW2PK7#h%>#VE#d9x}#UMsgDR9bU1ZRVE5OeNS@TnU;*z^ z-^p8aEVW{s4IqDCIUg=<$vvOMm~(bX^C!-#lIEMvS1(*y zxqqxPu+yt`2ybW15*qkN7`hqj{z4SZSlD}36dj=5pgfTFUQzUb+P)M;{(+@GFN*D; z4$#A(9?%}pUeL>2e^C@QxOa`9BS5`hiQ+_1y-yTpak)feqSy-R0qqH- zf~@2O^EpvWhlUjeOCLaXu>32l{S203v^gu!*b7k^c%QIm@P8tm&dP6?gUP9~iqG!F z|N3FrAh1?o3qkonBW4;T<1(;Opv)CvQ+*f`p(zfPKDX=~L&6>f;e*>x1U4TSr?AgZ zorQf#=yXXr3t-~~Y*bNRhO#u2*^4UFfE4~~VQL@D2h z;u8|W82=!&0owuW0>@OxB5YL=8W&hz|K*oO@c|X1KaF>)H$iq&-}sy$FSwIeDLFeW zXb@~~Kh4YjP0%8g)qR7mW?MpRnBD~ObH7*w@h=dU4pPO00?i=UR-6O6?uf#NT~iQw zwJWvq$0^g3GuJAtqy@QF!_hIAiGdffG#s=YIu~{1rVD%M9C!-}dFckX{H=&!q~RuK zp4E_nE6oT0isi#YEAz?{`{f@!T zhbi-*Je^NfMgY??#jS#I8o)xZXrj{?5Xv&B+*vH>sRNRZ*+5MehWE$t5x@-W=otFa zIW2Z#)^H_vX~dCm7MPF)XMr)mgp#a;A=QF@k=c||*fbgCE{y>w+2 zlraMiV;~sg3NgFYg47OR0(yKIZdF&dPRP z{OsdHNlu@Pdt_W{$k5-EaecC6&U%m#KGd|k#1xVk=wHZy{K-@!sUxQ|2$*0Koe))>74dZozU4;RwD?#6kDEd*rWElR;Y#M z0#*8wLQ`mGd3uP_mHLC|MEL}*M{b24GF0T;?T<>Ibkp3522H=;JEIL1Ci(skTaUX^pVzq|44pO1FhStSZJv~zZjwyngEWz*@|IU%Ryzie0p I`m`+Je^=5RI{*Lx diff --git a/examples/06_gamma_beta/example.c b/examples/06_gamma_beta/example.c index 3a7d422..00330f2 100644 --- a/examples/06_gamma_beta/example.c +++ b/examples/06_gamma_beta/example.c @@ -43,10 +43,3 @@ int main() free(seed); } -/* -Aggregation mechanisms: -- Quantiles (requires a sort) -- Sum -- Average -- Std -*/ diff --git a/examples/07_ci_beta/example b/examples/07_ci_beta/example new file mode 100755 index 0000000000000000000000000000000000000000..d5a34625d8c245d9bc6e48486d449e4b7d43aa40 GIT binary patch literal 22320 zcmeHPeQ;aVm4CLZ*hy?X1(ISuQbbS&^AXw3SL4z&Qer31V8A2}b%KMU*i!6`Z5dfg z9LlB!+et()n{L@9-O`TBbf=l!*>3tZleXI;J0CbesS{uuO4&F8%2S+xQ3wPG(f-c; zke_4;b{VFBbbQCsx%b?2&OP^>^Ul5RKI?s}qH1||R+b`*N%^`$tG60H%Mn39tE@s_|^NwUVlwudjM`+@!7JSu$ zPeMiTNhs|}HlfcVVb4@(Lmw3z`&N`xzFNVT3d;q(gfs?~Rc{6xKF@g7g0HZY%SnY7 zIa-lWj(0csgsBC^rvvHcm14Y1euxNrN=W{uQrz6MVd4DZ=7u@VO|7AhIUS1@&RIC$ z(e8K5W&NZVjX^QB;?C9NjtJp|8JRDpi#GKyUjNwAhxVJd7T)l$n;Kf?UH|wk-#Q~@ z(3nI+h4`h5sZI6E@gW*M|A$cHID$vLFSHQMRDOgGaVaN}{3!!|D)6=>O5XT6Ksq^% z8Ssy0z}IBJ*Jr>FX29P7ybYh!QVT>nJ-?p;|3C(Oa|Zm54ETk>+we&(_W_YkpXv0DPS^2{Z~rqcuSv>?Xop!lM|mG@D#YIgNVHU@H#HT|)lW+&;kZ0gLSFDn zlAoqnvD1Jz=K0+QJdGjC z9s^$5i2C{rc*^0!?B@%F-0S|{J%6g4yvQt?OXk0j5KI!CatY@!9hL z)N2)@Z={}B^kh69-_6>zA&5nfCfc+Sh(-GoZQ7{BqCZWvX+sc;{xH#|jX*5=?L?b4 z0I}#pi8d|#vFQDYHZAzEs3*~;4Nok3ccM)jomkYBXwya@7A;M*X#)_8UYBUo!XJxH zlWolQwXf4`%M~_Xv=^te7o@c3q_k(Jv}dNYXQZ?XQrhN}_9r+FFzRzYrF}M~{box0 zwUqYD#`fx`$)Jb=tKJJ6XuXT?fk(Ax6djJ!!iTl+A^jIupxV3O<(Z1|G`ZWSzH`K< zsoVN}>b7!w{(f>fa1pIPYCET{1^x-_iRo0d^gTCaH}fwBw183e z{g|fS6p5*8%w9|jBAl?)K-X{K4^#fv?_ZZ0!H0Q=M+ZG2lTZU<8hk8 zzk!%>t7{`_(Q(x=cuu|h)pP2!I(yER%@9{^&m+1=K&OYmL9>lE4Hh;}rpO+HFG9f2K=MC-}P37=-PSJVwU{k4zd@sTcr$ta#+D~v7tfgs3W zuiTDV#+jP?uOe&I!5VA&p~=9f`)lO+ME-$)iZCXuSZ=pMVP?)k9)KCf{!Cvzp6aXX znRd)q`F7^2Z;x};e}2SW6~8j0tMa#e6Ryg!lc$D3rauRsBv0k@`L8Lr7s3iOqppKR zum#Ugnf>++^4spr{#r}^nuF1z-$YHyWB&SU*za|uP|S>CJ)`~nA9DLW4cZK$ZnR$? z$n;&{`E2_|9|sw~{xik>|H*!5FQl>Gsu(UwE&MU!(8e!i{g;J64zW;I zr%i|F*Y@lGgnn(82_)0)=mqxxCh;uAESx@w%K%Y)pZNt^a7u?*WYfi6@Jm9Z$#~qV zhOC<|1zoQ0{C;X!S3#g~OBF(mQqXtD#q}RQ;C{IQ#2Nf@jk@Do;A3gX3V~)5X!W>#q1@4~`EL-N` z8Cc#P?!*=Wb}k5fiK%TSh{?UmHAF%#@~MvF#2QoYKBDdzW%SmWPk|f0c?Bk~vlm8x zgjr*KW$xu`$r|$qj0D+;=X5>ApyofKUw`<0Zq#FFfcq>3ky{fM0 z8z0sD3nr56rhgZVU`RGB*&K#N#35^}fy%W7AL6_p)BR9=cn1_fw6VaURrzxz#mm+3 z6)a$_D>{xkIq5|}TqkRA*;r@9aZ-r7`NyDvi*fB0AKh{vsGumEX45B~C)?9_mM18fDC&}omS*Ea z)_pEzVQqvE{oFZ{JEk23rQZmRJg?`%8qw`c5!#sOE4Xz~jebNMi?QtiyWd)^g+HO) zIxI@mVYXqR+=pE@a5;9eydL1_>2%Fo*I}-d;ozF00ppg7)}x;Hw^est$_iG8V_NtE zU0E<)|EU};>!aN;4a^&1wnAOIKt|Kb210)ck9G}T>iSC-Ho|*=%z1|GAMn+Q&|%l- z0ot{`6B#*Sn1I;-(&k}O*GE~~W_6v;-cL6aTGv1}u4^WVCUF>}X3;=rjq0-bL*H7` zx8eVr!p)KVrK6M6zWW5s!H(gG$t<#^O3lxEQ43vMuIIkh=zK#KKLEn9BOBrgk z5md-^MSmW7jx>pG#-%gn*Ed+tL7FO~Iv8D!e%iF71J^Gl5WV|dc3Gi{@S_fritCf? z=lX{0p33kGuGN*{Ph2&w@Jp*#Y2kAupy!zS*l9!qB17$d8x=C-JFl{Pwo^GIUdxJo zhCR!}!wAgsvSGD*2_$ff9=N0aA~ZvusqBA?l3b3xDje5;Fv_&N9EQQ#7(N6D63iJE znr5xsH7RlViyN!t`y^M_#jI_A@9M3wUt%}KJcgNHASwT22pWW$k$tDte1Zsn1_Aw1 z=#9$oBrNgPuCVSlabaF(m=y=^x@Q2BY>@8+CGc+mRNzxnm5 z@CUAVjaeC)wi~j7KLb5b%i=YQFUS<2uX=|~jV<|w_pWo*xz@SvbGc=lC(_ZlJn5%( zBXK{H-fliF{PWr^dISla-xra4$*<9S8S0%^Sx6m-R-h^GZ(TK&;UlZGhtA|-acah+P5DMk^g$jacw)m5yM4V~%}hvF}?MKaOPHz-!{tJ|%Zw6puH+MSCE z7*0V*!HZDk*JK2Q-s{K<&m$i^r+zP@Rt$7U)SiDq4=_j7?+vKs3`2L-rO1D@8-oy9 z;@RSf{ZkR@7>RVO!<^vTv89z!7hh0!{E0~{Mjql?oo4N#Q;>x{UhPRz>g&utOc#s{ zT@H0eLYG00R}y-3KY%zTPlBC!l0C(YxmfzqcEN}VH>y291r73$ll+KF;lNmYW3(V? zTKw+jI*g83ES@on^TN%Sp=ZI{r~z|M-4R4Qy7o_^DC#Yk4Rq zyN39SyB@TcKoaWsuCZHtpw{JV6ml(yk^P7m*&M_{do?2?Vsv8_?70NHpjBtD#dkAp ziE%oz2$p+I-9EtlTdb?w7xA!zhifVH%(a++2)#8}P8@?{VK;-3I))PN?lG|;D9E@j zRqTp@vTW+4R~6F>aR=E95AAF{SP7Z6djbX>vEfAQ(w4sXV|53=UnoAIZl|aHnC=ih zOv~tMauL}Ll=|ie$cv~IynqvD9wE=8224!d@f2L9_14+jC=WH>iP0z?STx3DDJ2|CQdaYB8>M~I&F_!e|z38Ed2-p*(d{c?t1#SEh}rG3Wm zwCG07Y<)3AI$K=~0?vuV4g7y=iWtrR!vZ)xE{t$!G|QT?o@V-4oD zVVv}#Hl^b>_>!GHA{nhQNAb8oKZxBgWdTI5Ki^&%Mh&R}Cn3NUCX5ui8Ck&w>rhT$ ze#Q`8C&hWCu1y{kCz`d9kss+L(3@S4U#z}mr0+{OFZrVN?R-nY{x|_#Lh5K%Oi960WOd7plreLT7%BP(=uR@Vmb4S(N^2 z5wN-m%)pW5@i~ITyLP+tNpdKbqam?7^IB} zHY2D;T=!vbYttqV;lV;Homk_#Due_{$=Tl%WFZmPnus_@(8;kR7z=_fEV?7&_NoyPEahCPvC#isG-$rnMA zNtFCa7)9$unXtfw1tu);w_|}U@n3soC9`MS>iw<0rUq|oz1N2S`1A%gdzx*2pKXIT z=qW0_W$ul(rgq!CvwgSNX8YD*?`9}YULRfJyalYb+3mKX*$p=;4o5Nm+V2lFZQR)G zb<|_0F6CmZApZTcsdb|**yy$4A67SdgEov9^4ePcL7UGXYHeV`e4*C*V3WTUf;W2u z?Ou0%gHI`(Q(p7&b(7_|w+dY%9Dk+24=?%6-))?{^6n|l)vUaXe`26A%h_`1($Bsb z>l5wKr-2-f$N%t*fGe`DHr;v=jW6SqcOo7y!S3Y(U@c%CDBA#az$XD4PsZbgfVHRM z@tuGjfDypGfG+_a2h;)cUXRBwQa_*-yRyA+#N)RB7Gl|Q0@h-YtOMMGrImJVi?Gf< zN$oT7`161|;2FT$v++3n^XDGG>41f}=_~@Q#a($h;7**q>;}|g;G;H<6Y)>5N+S;I z+W-skFywK-LBRchJF)k;0C)k=DgG%F=dP?|J+M;A>ab>Aou8Yx8@dx7Pi++qb1x=} z@?3V~Wk+>BK0~j?T=mtP=GkWvn)vJR+4DcJ0U<;od=Q`A zzz+)?0it;npFxa=$9JsQ@ZE<`IsCB&P^qw3w`W&cY}sCmb(w{BD~#_T;6CtOM+g*> zbsV4PfxD3qD1QJ{^yfU#SqT5V!|{+sZI}o)p76aJQ1arkPn5+ZyHs4NVPw`G{Mh zQI=b*|B$`ZQuvsu!eZ;qsjw71YA&~wY@g(^Eb7W#VW~FNVv>tomJ*kxXsN{p(xn#b zQcIpTbrIMfhm29g{LfJ*S&wFyTMD~OD=d9E_gVJl+-BLA699Cfy>^PFu$-~N03{%& z-`fOnozQib;N5O=S!`W7D=h0x4`N^&#JL#L?NggUmIK>%;)L+H;OQmZA2pR*Y}+wN zQJ48POP|@)mu1;&uCVMgHv?j9bN0F^mLd?jr1q|<*_BfS(QQ-HRV$|2FkTSjjpCZ! zhdS4B1%hS^99hBQjA-_N=0QAK$VBrs(EJx@@^H=n26dA280p<>fm(P{%(u5KvYs;8%Hbo0${={kX+pZ8;Av{fqlf{7ySXWQmm{oB~!0xL&|E0XqfUEnuI3 z`vi;#I3(b(fTIE`#vh-h=9?7S(&9~u{QbF6I{miu!L;<|M4TtnPf}!jC*yM!c^xF< z^As5$$@s~NjH_h)6h+2wGJdKe<2)IUt6uuY&|$*WZ%d1}DAqK7G%2cLON)n)G(Rph zAsA#{F${~N0}z)OF}!Y4E?4C5@r}JE+_e|ZHDVxY zPidR*GXJB=iVNRa20cZ3WYft?V;?%lAUSe9;ZpFQFig&PJ3j?{I{h=v>pz&B=}PlL zDR~U0lf&*Ek^SX5BdMQfa>o1dWkxVwPo3e%$JNJ7&Un9{N8-+=RT3d>EYo1X%Yd|! z!BOP8b{)fWobh{c2K-$@PPJG_+4%_=>wvcz1&I14=wnmH$6*JPGoE1=!&{XSahSAE zko-X4cjJEpsa(f_@)O{#ke?dQzhv}+$N;V;eD?V3`cMY?lMIh|77HpnQ-Zt^fiIze zM`7nsfFBC{`s*25af*dBA222K8NoFSZ&h}Sg`b^`p{F$c(w{g;>?;hZ%WnZbt^-OmDBjIV}5{|*q1S>xTp5m%aY&7Iq*H0wg#_8Dk z6{^1#_|77Zm;U{|z>f<5%Ckf|J4!d-5kcQ~6PM%VK>0WWy^8gZ=_xi;(qD97ZAG3b ztQGmb6fDVF}z?lfJcBQf9)6klKa2j4D>(9fZs=Q5D%&P&>0!&9}yeeO=9fl zfVUa+kkvOc$p1jdVGjyGEd7{;0Au!<$t6?>K?MTeH=85a*&?blffp_jlpK^X`Xm%< zR0=t5!Y{Z&SXnFZHgP>}5O}<~!_(5np25@yg6+W&o_W_RV{Z|02V30r^d6CRyrQGQ z@7~z#-{5I>Hw68GcDEt?1OlF|Zf|QauoW-&@U(c{ z4WX8ntza><+@K65ks6*7x!ucGy6&iOSKL`nFEL4`ZBX3hYwmR2QMoL+hrROzD7;C> zU7-mHTKP)Fy`t)_rLHRXUCWoRs;F_-xRzE`kW}^}lKS=#(;98#%S)V2`9>0Y%Lu)q zWDG;^C~0TLbZ_*uw0J)Cg)HnPDCu|u^cV}44$*5{$he*W48@vofw#dE^w7.J3O zW73OP#yaOqje>qzYxM{4UKYjh&`geo7rl&e&MlJ+A@YdeqIbaFia)l}XEs)nXkcc|UlASBd>TG~S`Qlplp zj$kO@6-{oj4gL6` z$De>IZ2^D1x4qr%4Fvpw#MtxZ3YA;Do>okB11NY1vbWl}4R6gE8=8Wtx{_Y7GlrvZ zqSxvqW87}M*GF6cLRt972cLe*OMS=~$zp_G#!O1ja+6|$?3w(Io@5DlQIG*qrjP)4 zse7)QYZZ{bV$ma2R2=PFTY{bqfWZKVjk3n<+NL=0kk{+j*cx)czHQz>aH|2eA=HF7 zZ#6Xtlxu0_oSxz4eMC=?(A48%5@7Y=AV`;xc#&;~cV%p*$^35XRXN8Hxi} zcneaZ%*LM^(j5-AXxJn3BJZyyl;=UBo6*WT!Rzs%qZi39@53duq0d<4x+`&X`(Wg+ z#rj7@LOX9JzOfu7l}faX{PMnDLR%6w-Qy&eZ=i|aW6^en6?x7jp&VbHbEWF<0T7lm}Bs^pjFQBA4heqge|bM7VIc`O6zM;S z`aLQ{BKhTck%Si|p@6C5znaQFmlMkKrP>sR)bZa0IlYfdq@Be1lYGBpD!pt=_%3)= z1HU}q&=v!UBI$T*tS#}Mpo3ypjxW!% +#include +#include + +// Estimate functions +double beta_1_2_sampler(uint64_t* seed){ + return sample_beta(1, 2.0, seed); +} + +int main() +{ + // set randomness seed + uint64_t* seed = malloc(sizeof(uint64_t)); + *seed = 1000; // xorshift can't start with 0 + + struct c_i beta_1_2_ci_90 = get_90_confidence_interval(beta_1_2_sampler, seed); + printf("90%% confidence interval of beta(1,2) is [%f, %f]\n", beta_1_2_ci_90.low, beta_1_2_ci_90.high); + + free(seed); +} diff --git a/examples/07_ci_beta/makefile b/examples/07_ci_beta/makefile new file mode 100644 index 0000000..ef385f7 --- /dev/null +++ b/examples/07_ci_beta/makefile @@ -0,0 +1,53 @@ +# Interface: +# make +# make build +# make format +# make run + +# Compiler +CC=gcc +# CC=tcc # <= faster compilation + +# Main file +SRC=example.c ../../squiggle.c +OUTPUT=example + +## Dependencies +MATH=-lm + +## Flags +DEBUG= #'-g' +STANDARD=-std=c99 +WARNINGS=-Wall +OPTIMIZED=-O3 #-Ofast +# OPENMP=-fopenmp + +## Formatter +STYLE_BLUEPRINT=webkit +FORMATTER=clang-format -i -style=$(STYLE_BLUEPRINT) + +## make build +build: $(SRC) + $(CC) $(OPTIMIZED) $(DEBUG) $(SRC) $(MATH) -o $(OUTPUT) + +format: $(SRC) + $(FORMATTER) $(SRC) + +run: $(SRC) $(OUTPUT) + OMP_NUM_THREADS=1 ./$(OUTPUT) && echo + +time-linux: + @echo "Requires /bin/time, found on GNU/Linux systems" && echo + + @echo "Running 100x and taking avg time $(OUTPUT)" + @t=$$(/usr/bin/time -f "%e" -p bash -c 'for i in {1..100}; do $(OUTPUT); done' 2>&1 >/dev/null | grep real | awk '{print $$2}' ); echo "scale=2; 1000 * $$t / 100" | bc | sed "s|^|Time using 1 thread: |" | sed 's|$$|ms|' && echo + +## Profiling + +profile-linux: + echo "Requires perf, which depends on the kernel version, and might be in linux-tools package or similar" + echo "Must be run as sudo" + $(CC) $(SRC) $(MATH) -o $(OUTPUT) + sudo perf record ./$(OUTPUT) + sudo perf report + rm perf.data diff --git a/makefile b/makefile index 1a81168..5664c73 100644 --- a/makefile +++ b/makefile @@ -11,6 +11,7 @@ all: cd examples/04_sample_from_cdf_simple && make && echo cd examples/05_sample_from_cdf_beta && make && echo cd examples/06_gamma_beta && make && echo + cd examples/07_ci_beta && make && echo format: squiggle.c squiggle.h $(FORMATTER) squiggle.c diff --git a/squiggle.c b/squiggle.c index 96288d6..3ddb416 100644 --- a/squiggle.c +++ b/squiggle.c @@ -94,6 +94,12 @@ double sample_gamma(double alpha, uint64_t* seed) // A Simple Method for Generating Gamma Variables, Marsaglia and Wan Tsang, 2001 // https://dl.acm.org/doi/pdf/10.1145/358407.358414 // see also the references/ folder + // Note that the Wikipedia page for the gamma distribution includes a scaling parameter + // k or beta + // https://en.wikipedia.org/wiki/Gamma_distribution + // such that gamma_k(alpha, k) = k * gamma(alpha) + // or gamma_beta(alpha, beta) = gamma(alpha) / beta + // So far I have not needed to use this, and thus the second parameter is by default 1. if (alpha >= 1) { double d, c, x, v, u; d = alpha - 1.0 / 3.0; @@ -377,6 +383,43 @@ struct box sampler_cdf_double(double cdf(double), uint64_t* seed) return result; } +// Get confidence intervals, given a sampler + +struct c_i { + float low; + float high; +}; +int compare_doubles(const void *p, const void *q) { + // https://wikiless.esmailelbob.xyz/wiki/Qsort?lang=en + double x = *(const double *)p; + double y = *(const double *)q; + + /* Avoid return x - y, which can cause undefined behaviour + because of signed integer overflow. */ + if (x < y) + return -1; // Return -1 if you want ascending, 1 if you want descending order. + else if (x > y) + return 1; // Return 1 if you want ascending, -1 if you want descending order. + + return 0; +} +struct c_i get_90_confidence_interval(double (*sampler)(uint64_t*), uint64_t* seed){ + int n = 100 * 1000; + double* samples_array = malloc(n * sizeof(double)); + for(int i=0; i