Dün konuştuğumuz Monte Carlo tümlev alma (http://ddili.org/forum/post/10353) yöntemi hakkında birisi bir Python programı yazmış. Sonra Riemann yönteminin daha hızlı ve daha doğru sonuçlar verdiğinden bahsetmiş.

Ben de hızlıca bir D kodu yazayım dedim.

import std.stdio, std.random;

void main()
{
   int icerdekiNoktalar = 0, toplamNokta = 0;

   Random cekirdek;
   for (int i = 0; i < 100; ++i)
   {
       auto x = uniform(0.0L, 1.0L, cekirdek);
       auto y = uniform(0.0L, 1.0L, cekirdek);

       ++toplamNokta;
       writeln("toplam nokta", toplamNokta);

       writeln("x,y ", x," ", y);


       if (((x * x) + (y * y)) <= 1)
       {
           icerdekiNoktalar++;
           writeln(icerdekiNoktalar);
       }
   }

   writeln("icerdeki ", icerdekiNoktalar);
   writeln("toplam ", toplamNokta);

   float pisayisi = cast (float) (icerdekiNoktalar/toplamNokta);


   writeln("Pi sayisi = ", pisayisi);

}

Ama nedense pi sayısının çeyreğini bulurken sonuç 0 çıkıyor.

Sanırım biraz da D çalışmanın zamanı gelmiş ;)

--
[ Bu gönderi, http://ddili.org/forum'dan dönüştürülmüştür. ]

Sanırım şu şekilde yapmalısın:

   float pisayisi = cast(float)icerdekiNoktalar/cast(float)toplamNokta;

Ama belki de hiç dönüştürmeye gerek yoktur. Bütün değişkenler double olamaz mı?

--
[ Bu gönderi, http://ddili.org/forum'dan dönüştürülmüştür. ]

Oluyormuş:

/*
* mcPi.d: Monte Carlo Pi
* Source: http://ddili.org/forum/thread/1292
* Programmer: Erdem Öncel
*/
import std.stdio, std.random;

void main()
{
   double icerdekiNoktalar = 0.0, toplamNokta = 0.0;
   double resolution = 2000.0;

   Random cekirdek;
   do {
       auto x = uniform(0.0L, 1.0L, cekirdek);
       auto y = uniform(0.0L, 1.0L, cekirdek);

       toplamNokta++;
       debug(1)
       {
           writeln("toplam nokta", toplamNokta);
           writeln("x, y ", x," ", y);
       }

       if (((x * x) + (y * y)) <= 1)
       {
           icerdekiNoktalar++;
           debug(2) writeln(icerdekiNoktalar);
       }
   } while(icerdekiNoktalar < resolution);

   writefln("icerdeki: %.0f", icerdekiNoktalar);
   writefln("toplam: %.0f", toplamNokta);

   auto pisayisi = (icerdekiNoktalar/toplamNokta) * 4.0;

   writeln("Pi sayisi = ", pisayisi);
}

--
[ Bu gönderi, http://ddili.org/forum'dan dönüştürülmüştür. ]

Teşekkürler double kullanınca çalıştı.

Ama üçüncü basamağa kadar doğru tahmin edebilmesi için bir milyondan fazla tahmin yapması gerekiyor.

Hatta aslında sorduğum soru buydu (http://temelelektronik.net/sorular/115/%24-displaystyle-int_0-1-sqrt-16-16x-2-%24). Birisi videoyla cevap vermiş.

http://www.youtube.com/watch?v=djJzd-MtVgM

Orada ilk anlattığı benim çözümüm ;-) Ama ikinci çözüm de güzel.

--
[ Bu gönderi, http://ddili.org/forum'dan dönüştürülmüştür. ]

Bende bir kaç milyonu denedim ve hislerim Pi'den uzaklaştığını söylüyor...

Belki de beşik gibi bir salınımı olabilir. Çünkü bir ara neredeyse 3.14'ün altına indiğini ve sonra yükseldiğini gördüm. Belki belli aralıklarla numune alınıp (örneğin 100 bin'de 1) değerler bir çizelgeye aktarılabilir. Çıkan eğri çok ilginç olabilir...

--
[ Bu gönderi, http://ddili.org/forum'dan dönüştürülmüştür. ]

Çok güzel program. Teşekkürler! :)

Ali beyin yazdığı programı 10 milyon sayıyla denedim üçüncü ondalık basamağa kadar buldu.

Python programını yazan kişininkini de denedim. O da üçüncü ondalık basamağa kadar buldu.

Pi sayısı 3.14159 şeklinde devam ediyor. Gerçekten de Riemann yöntemi daha hızlı çalışıyor ama daha hassas bir sonuç bulamıyor.

Merak edenler için Python programını da yazayım.

import math
import random

# from time import perf_counter
# burayı alttaki satırla değiştirdim

from timeit import default_timer as perf_counter

def monte_carlo(start, stop, num_points, f):
   hits = 0
   upp_bound = 0
   for i in range(num_points):
       x = random.uniform(start, stop)
       y = random.uniform(0, 4) # this is cheating since I already know it's 4
       if y <= f(x):
           hits += 1
       if y > upp_bound:
           upp_bound = y
   ans = (hits / num_points) * ((stop - start) * upp_bound)
   return ans

def riemann(start, stop, step, f):
   return sum(f(start+step*m)*step for m in range(int((stop - start)/step)))

f = lambda x: math.sqrt(16 - 16*x**2)

start = perf_counter()
a = monte_carlo(0, 1, 10000000, f)
m_time = perf_counter() - start
start = perf_counter()
b = riemann(0, 1, 1/15000, f)
r_time = perf_counter() - start
print('Monte Carlo: %s\nTime elapsed: %s' % (a, m_time))
print('Riemann sum: %s\nTime elapsed: %s' % (b, r_time))

Burada işleve bir fonksiyon göndererek direkt olarak bahsettiğim sorunun integralini hesaplatmış anladığım kadarıyla.

--
[ Bu gönderi, http://ddili.org/forum'dan dönüştürülmüştür. ]

Bende istatistiksel açıdan inceleyim dedim...:)

Malum hissettiğim bir düzensizlik vardı. Bunun matematiksel karşılığı nedir bilmiyorum ama bazen +/- dalgalanma (peak) yapıyor. Tabi sonuçları, örnekleme çözünürlüğü on binde bir olacak şekilde çizelgeye yansıttığım için derin/yüksek tepeler pek belli olmuyor. Kullandığım kod ise aşağıda ve bunu döngü içindeki koşula (denklem kapsamına) koymalısınız...

           enum çözünürlük = 10_000.0;
           if(içerdekiNoktalar % çözünürlük == 0) {
             auto örnekleme = (3.141592 * çözünürlük) -
                              (içerdekiNoktalar/toplamNokta * 4.0 * çözünürlük);
             //if(örnekleme > - çözünürlük && örnekleme < çözünürlük)
             // arada dalgalanma (peak) var!/neden düzensiz olabilir ki?
             writefln("%.0f", örnekleme);

           }

'http://img27.imageshack.us/img27/1205/lcv3.png'

Bu çizelgeyi, elimden geldiğince yumuşatmaya ve gereksiz ayrıntılardan uzak tutmaya çalıştım. Gri kutu ile işaretlediğim başlangıçtaki kararsızlık normal çünkü henüz iç/dış noktaların oranı yeterince büyük sayı değil. Ama sonradan 1 milyon döngüye rağmen, bazen Pi sayısına yaklaşıyor, bazense uzaklaşıyor!

Dip Not: Tam olarak denk geldiği noktayı kestirmek güç olduğu için başka algoritmalar ile desteklemeden Pi sayısını bulmak için kullanılabilir gözükmüyor. Ama rasgeleliğin sonuçlarını değerlendirmek açısından güzel...:)

Sevgiler, saygılar...

--
[ Bu gönderi, http://ddili.org/forum'dan dönüştürülmüştür. ]

Bir tane de yalnızca Phobos ve aralıklarla:

import std.stdio;
import std.range;
import std.algorithm;
import std.typecons;
import std.random;
import std.conv;

void main()
{
   enum size_t toplam = 1_000_000;

   const size_t içerdekiler =
       iota(toplam)
       .map!(_ => tuple(uniform(0.0L, 1.0L), uniform(0.0L, 1.0L)))
       .count!(a => ((a[0] * a[0]) + (a[1] * a[1])) < 1);

   const pi = 4.0 * içerdekiler / toplam;

   /* Başka yollarla da olur. Bir de std.conv.to ile:
    *
    *   const pi = içerdekiler.to!double / toplam * 4;
    */

   writeln(pi);
}

Ali

--
[ Bu gönderi, http://ddili.org/forum'dan dönüştürülmüştür. ]