Buffon needle (Java and MATLAB)

 

 

Java code

 

/*************************************** * ************************** Java Applet buff. java March 1998-Michael J. hurben this applet provides a _ simulation _ of the famous Buffon's needle experiment. buffon's needle is a Monte Carlo Simulation for the determination of the value of pi. **************************************** * *************************/import Java. applet. *; import Java. AWT. *; Public Class buff extends applet implements runnable {dimension offdimension, D; // variables used to create an image offimage; // offscreen image via the update () graphics offgraphics; // method, to reduce flicker int needlelength = 20; int numrows = 10; int xborder = 15; // basic layout variables int yborder = 45; int margin = 50; int tableheight = needlelength * numrows; // int tablewidth = 200; // def INE size of table int tablexcorner = xborder; // where needles are thrown int tableycorner = yborder; // int barwidth = 20; // define size of bar graph int barheight = tableheight; int barxcorner = xborder + tablexcorner + tablewidth + margin; int barycorner = tableycorner; double bartop = 3.24; // set limits for limits of bar graph double Barbot = 3.04; double barmid = (bartop + Barbot)/2; int barvalue; // Dynamic value of the bar graph int plotwidth = 200; // define size of plot int plotheight = tableheight; int plotxcorner = barxcorner + margin + barwidth + xborder; int plotycorner = tableycorner; double plotmax = 4; // set limits on the plot double plotmin = 2; double plotmid = (plotmax + plotmin)/2; int barmaxplot = (INT) (plotheight * (plotmax-bartop)/(plotmax-plotmin); int barminplot = (INT) (ploth Eight * (plotmax-Barbot)/(plotmax-plotmin); int hire ick = (barmaxplot = barminplot? 1: barminplot-barmaxplot); int pointerwidth = 30; // bar graph pointer int pointertop, pointerbot; int power = 5; // orders of magnqueue for the x-axis of the plot int n = 0; // The number of needles thrown int hit = 0; // The number of needles which hit lines int I = 1; // double needlex1, needlex2, needley1, needley2; double deltax, deltay; int ysign; final double pivalue = 2.0 * Math. atan( 1.0); int X1, x2, Y1, Y2; // used to draw needles int GX1, GX2, gy1, gy2; // used to make the plot int speed = 20; // sets speed with try-Catch Block double H = 0; // double version of the number of hits double estpi, prevestpi; // estimates of PI double prob = 1; // ratio of hits to trys double Prev = 1; // ratio on the previous trial Boolean oldscreen = false; Boolean begin = false; thread t; button B1, B 2, B3, B4; Public void Init () {setlayout (New flowlayout (flowlayout. left); b1 = new button ("START"); b2 = new button ("slower"); B3 = new button ("faster "); b4 = new button ("done"); add (B1); add (B2); add (B3); add (B4); t = new thread (this); T. start ();} public Boolean action (event E, object O) {If (O. equals ("START") {n = 0; hit = 0; oldscreen = false; begin = true;} else if (O. equals ("done ")){ T. Stop ();} else if (O. Equals ("faster") {speed = (speed> 10? Speed-10: 10);} else if (O. equals ("slower") {speed = speed + 10;} return true;} public void run () {While (true) {If (BEGIN) {// calculate the position of the needle ends... needlex1 = tablewidth * Math. random (); needley1 = tableheight * Math. random (); deltax = needlelength * Math. sin (2 * pivalue * Math. random ()-pivalue); deltay = math. SQRT (needlelength * needlelength-deltax * deltax); needlex2 = needlex1 + d Eltax; ysign = (math. Random () <0.5? -1: 1); needley2 = needley1 + ysign * deltay; // check to see if it crosses a line... for (INT yline = 0; yline <= needlelength * numrows; yline + = needlelength) {If (needley1 <= yline & needley2> = yline) | (needley1> = yline & needley2 <= yline) {hit ++; break ;}} n ++; Prev = prob; H = hit; prob = H/N; prevestpi = 2/Prev; estpi = 2/prob; // calculation for the plot GX1 = (INT) (plotxcorner + math. log (n) * plotwid Th/(Power * Math. log (10); GX2 = (INT) (plotxcorner + math. log (n + 1) * plotwidth/(Power * Math. log (10); gy1 = (INT) (plotycorner + (plotheight/(plotmin-plotmax) * (prevestpi-plotmax); If (gy1 <= plotycorner) gy1 = plotycorner + 1; if (gy1 = plotycorner + plotheight) gy1 = gy1-1; gy2 = (INT) (plotycorner + (plotheight/(plotmin-plotmax )) * (estpi-plotmax); If (gy2 <= plotycorner) gy2 = plotycorner + 1; if (gy2 = plotycorner + Plotheight) gy2 = gy2-1; // Calculation for drawing the needles X1 = (INT) (xborder + needlex1); x2 = (INT) (xborder + needlex2); Y1 = (INT) (yborder + needley1); y2 = (INT) (yborder + needley2); repaint () ;}try {thread. currentthread (). sleep (speed);} catch (interruptedexception e) {}} public void paint (Graphics g) {d = size (); Update (g );} public void Update (Graphics g) {If (offgraphics = NULL) // se Tup an off-screen image | (D. Width! = Offdimension. width) // via the update () method. | (D. height! = Offdimension. Height) // {offdimension = D; offimage = createimage (D. Width, D. Height); offgraphics = offimage. getgraphics ();} If (! Oldscreen) {offgraphics. setcolor (getbackground (); offgraphics. fillrect (0, 0, D. width, D. height); offgraphics. setcolor (color. gray); offgraphics. fillrect (plotxcorner, plotycorner + barmaxplot, plotwidth, barminplot-barmaxplot); offgraphics. drawstring ("+ barmid +"-", barXcorner-45, barycorner + barheight/2); offgraphics. drawstring ("+ bartop +"-", barXcorner-45, barycorner + 5); offgraphics. drawstring ("" + B Arbot + "-", barXcorner-45, barycorner + 3 + barheight); offgraphics. drawstring ("-" + plotmax, plotxcorner + plotwidth + 2, plotycorner + 5); offgraphics. drawstring ("-" + plotmin, plotxcorner + plotwidth + 2, plotycorner + plotheight + 5); offgraphics. drawstring ("-" + plotmid, plotxcorner + plotwidth + 2, plotycorner + plotheight/2 + 5); offgraphics. drawstring ("estimate. log (number of tries) ", plotxcorner + 5, plotycorner + plothei Ght + 20); // draw lines which relate the bar graph and the plot... offgraphics. drawline (barxcorner + barwidth + 5, barycorner, plotXcorner-5, barmaxplot + plotycorner); offgraphics. drawline (barxcorner + barwidth + 5, barycorner + barheight, plotXcorner-5, barminplot + plotycorner); // draw tabletop and the plot axes offgraphics. setcolor (color. blue); offgraphics. drawrect (tablexcorner, tableycorner, tablewidth, Tab Leheight); offgraphics. drawrect (plotxcorner, plotycorner, plotwidth, plotheight); for (I = 1; I <tableheight/needlelength; I ++) {offgraphics. drawline (tablexcorner, tableycorner + needlelength * I, tablexcorner + tablewidth, tableycorner + needlelength * I);} oldscreen = true;} else {offgraphics. setcolor (getbackground (); offgraphics. fillrect (0, 0, D. width, yBorder-10);} // write N and estpi to screen offgrap Hics. setcolor (color. black); offgraphics. drawstring ("current estimate of Pi:" + estpi, plotxcorner + 10, 20); offgraphics. drawstring ("Number of tries:" + N, 220, 20); offgraphics. setcolor (color. black); // make sure that none of the needles appear outside // The table boundary: X1 = (x1> tablexcorner + tablewidth? Tablexcorner + tablewidth: X1); x2 = (X2> tablexcorner + tablewidth? Tablexcorner + tablewidth: X2); Y1 = (Y1> tableycorner + tableheight? Tableycorner + tableheight: Y1); y2 = (Y2> tableycorner + tableheight? Tableycorner + tableheight: Y2); X1 = (x1barheight + pointerwidth/2) barvalue = barheight + pointerwidth/2; pointertop = (barvalue-pointerwidth/2 <0? 0: barvalue-pointerwidth/2); pointerbot = (barvalue + pointerwidth/2> barheight? Barheight: barvalue + pointerwidth/2); offgraphics. fillrect (barxcorner, barycorner + pointertop, barwidth, pointerbot-pointertop); // center line in the pointer offgraphics. setcolor (color. gray); If (barvalue <0) barvalue = 0; If (barvalue> barheight) barvalue = barheight; offgraphics. drawline (barxcorner, barycorner + barvalue, barxcorner + barWidth-1, barycorner + barvalue); // Finally, draw the plot offgraphics. setcolor (color. red); offgraphics. drawline (GX1, gy1, GX2, gy2); G. drawimage (offimage, 0, 0, this );}}

MATLAB

% Buffon test % for get the PI % Aris 2006.3.25% L is the length of needel % N is the numtest function buffonpi = buffon_test (L, n) winnum = 0; for I = 1: N center_l = Rand * l/2; % the location of the center of the L sta = Rand * PI/2; temph = (l/2) * sin (STA); If temph> center_l winnum = winnum + 1; end buffonpi = 2 * n/winnum;