Почему игольчатый метод Буффона при оценке числа пи не дает ожидаемого значения?

Question

Вот мой код

buffon <- function(a,l)
    {
      n<-0
      N<-0
      repeat {N<-N+1
      print(N)
      p<-c(n/N)
      # Sample the location of the needle's centre.
      x<-runif(1,min = 0,max =a/2)
      print(x)
      # Sample angle of needle with respect to lines.
      theta<-runif(1, 0, pi/2)
      print(theta)
      # Does the needle cross a line?
      k<-l/2*sin(theta)
      ifelse(x<=k,n<-n+1,n<-n)
      p<-c(n/N)
      print(p)
      pie<-(2*l)/(p*a)
      print(pie)
      if(N>5000) {break}
      }
    }

Я пытаюсь оценить значение пи, используя идею иглы Буффона, однако, когда я пытаюсь: buffon(2,3), окончательная оценка составляет 3.8 , что намного больше 3,1 .Может ли кто-нибудь объяснить мне, есть ли какие-либо ошибки в моем коде, или я не могу использовать пи для оценки пи?

---

Добавление:

Я понял, что многие строки моего кода являются избыточными, поэтому янемного изменил его сегодня утром:

  buffon01 <- function(n,a,l)
  {
    # Sample the location of the needle's centre.
    x<-runif(n,min = 0,max =a/2)
    # Sample angle of needle with respect to lines.
    theta<-runif(n, 0, pi/2)
    # Does the needle cross a line?
    k<-l/2*sin(theta)
    # l is the length of the needle
    # a is the distance between to parallel line
    v<-length(x[x<=k])
    p<-c(v/n)
    pie<-(2*l)/(p*a)
    list("pi"=pie,"error"=abs(pie-pi))
  }

Установив путь больше l, я смогу получить довольно близкий результат к 3.14 ... но полученный результат очень нестабилен ... как весли я снова проведу тот же эксперимент, 3.1 * может быть любым числом, кроме 4 ... Я проигнорировал некоторые другие проблемы в моей настройке?

Answer 1

Насколько я могу судить.Ваша проблема в том, что вы используете «длинную» иглу (где l> a), и тогда формула l/2*cos(theta) не работает.В этом случае вам нужно использовать немного более сложную формулу .

Так что я немного почистил ваш код и изменил, убедившись, что l

buffon <- function(a,l) {
  stopifnot(l < a)  
  n<-0
  N<-0
  repeat {
    N <- N+1

    # Sample the location of the needle's centre.
    x <- runif(1, min=0, max=a/2)

    # Sample angle of needle with respect to lines.
    theta <- runif(1, min=0, max=pi/2)

    # Does the needle cross a line?
    n <- ifelse(x <= l/2*cos(theta), n + 1, n)

    # Estimate probability of crossing line
    p <- n/N

    # Compute pi
    pie <- (2*l)/(p*a)

    if (N>50000) { # Note the increased iterations
      break
    }
  }
  return(pie)
}

ans <- buffon(a=3,l=2)
print(ans)
#[1] 3.159621

Похоже, что вы также переключили cos с sin в указанной формуле для проверки, попала ли стрелка в линию.Однако - вне головы - это не имеет значения.Наконец, печать внутри цикла «repeat» заставляет функцию печатать каждый раз, когда она встречается (и, таким образом, заполняет консоль).Вместо этого я изменил функцию, чтобы она возвращала объект pie, чтобы вы могли сохранить его в переменной.

Надеюсь, это поможет.

Редактировать:

Ваша новая краткая функция иллюстрирует проблему.Сравните:

buffon01(1e6, a=3, l=2)
#$`pi`
#[1] 3.143023
#
#$error
#[1] 0.00143062

buffon01(1e6, a=2, l=3)
#`pi`
#[1] 3.850819

#$error
#[1] 0.7092266

Здесь вы видите, что использование l> a дает сбой, потому что используется неправильная формула.Что касается сходимости к пи, игла Буффона довольно медленно сходится, и для получения достойной оценки требуется много бросков.Но это статистический вопрос, который затрагивается здесь

.