Introdução


\(\pi\) é a constante matemática mais conhecida no mundo, tendo aplicações em diversas áreas do conhecimento.

Originada a partir da proporção numérica existente entre o perímetro de uma circunferência e seu respectivo diâmetro, teve sua representação adotada pela letra grega \(\pi\) (pi), por volta do século XVIII.

O valor de \(\pi\) pertence ao conjunto dos números irracionais, sendo assim, sua representação é infinita e não-periódica, de forma que nunca se terá o valor exato de \(\pi\), o máximo que se consegue são estimações de seu valor. Uma aproximação com 10 casas decimais de \(\pi\) é dada por \(3,141592653589793\).

Nesse contexto de aproximações e estimações do valor de \(\pi\), surge como uma opção viável o método de Monte Carlo (MMC). Este método, se baseia em amostragens aleatórias massivas, para se obter resultados númericos — probabilidades — heuristicamente.

Para estimação de seu valor, utiliza-se do fato que \(\pi\) satisfaz a seguinte relação matemática:

Seja \(\Gamma\) uma circunferência, cujo raio é \(r\).
Denotemos por \(A_{\Gamma}\) a área de \(\Gamma\).
Então, a seguinte relação matemática existe: \[A_{\Gamma} = \pi \cdot r^2\]
Através deste resultado, a simulação de Monte Carlo é elaborada da seguinte maneira:


Segundo uma distribuição uniforme no intervalo \([-1;1]\), gera-se a coordenada no eixo \(x\) do ponto.

Segundo uma distribuição uniforme no intervalo \([-1;1]\), gera-se a coordenada no eixo \(y\) do ponto.

Repete-se esse processo \(N\) vezes.


Desta forma, gera-se \(N\) pontos uniformemente distribuídos dentro do quadrado \([-1;1]\times[-1;1]\).

Seja \(N_c\) a quantidade de pontos que se localizam dentro da circunferência de raio \(1\) e centro em \((0,0)\).

Denota-se por \(A_q\) a área do quadrado \([-1;1]\times[-1;1]\) e por \(A_c\) a área da circunferência centrada em \((0,0)\) e de raio \(1\).

Tem-se, então, que: \[\begin{cases} A_q = (1-(-1))^2 = 2^2 = 4 \\ A_c = \pi\cdot r^2 = \pi\cdot1^2 = \pi \end{cases}\]

Logo, a razão de pontos dentro da circunferência em relação ao total de pontos deverá tender à razão entre a área da circunferência e a área do quadrado. Então:

\[ \dfrac{N_c}{N} \rightarrow \dfrac{A_c}{A_q} = \dfrac{\pi}{4} \]

De forma que:

\[ \dfrac{4\cdot N_c}{N} \rightarrow \pi \]

E, é com base nesta convergência, que se estimará o valor de \(\pi\).

Desenvolvimento


Foram elaboradas duas maneiras distintas de se realizar o método de Monte Carlo: uma com a utilização de funções de loop e outra através do uso de vetores. As duas formas serão exibidas a seguir:

Método de Monte Carlo - Loop


Gerou-se uma amostra de tamanho \(N = 1000\) e, a partir dos dados gerados, foi estimado o valor de \(\pi\):

# ET591-EstComp        
# MatheusLeite

# Número de iterações
N <- 1000
# Variável que acumula total de pontos dentro da circunferência
totalDentro <- 0
# Inicialização da semente
set.seed("2507")
# Loop
for (i in 1:N) {
  # Gera-se dois valores segundo uma distribuição uniforme no intervalo [-0.5;0.5] 
  # Eles caracterizarão as coordenadas no eixo x e y do ponto
  u.x <- runif(n = 1, min = -1, max = 1)
  u.y <- runif(n = 1, min = -1, max = 1)
  # Verifica se ponto está dentro da circunferência
  dentro <-  u.x^2 + u.y^2 <= 1
  # Adiciona 1 caso a verificação anterior seja verdadeira
  totalDentro <- totalDentro + dentro
}
# Apaga variável i
rm(i)
# Estimação do valor de pi
piHat <- 4*totalDentro/N
# Imprime-se valor estimado de pi na tela
cat("O valor estimado de pi foi:", round(piHat, digits = 5))
## O valor estimado de pi foi: 3.136

Método de Monte Carlo - Vetor


Gerou-se uma amostra de tamanho \(N = 1000\) e, a partir dos dados gerados, foi estimado o valor de \(\pi\):

# ET591-EstComp        
# MatheusLeite

# Biblioteca Necessária
library(tidyverse)
# Número de iterações
N <- 1000
# Inicialização da semente
set.seed("2507")
# Gera-se N pontos uniformemente distribuídos em [-1;1]x[-1;1]
# Para cada um deles é verificado se encontra-se dentro ou não da circunferência
base <- tibble(Order = 1:N,
               U.x = runif(N, -1, 1),
               U.y = runif(N, -1, 1),
               Condicao = U.x^2 + U.y^2 <= 1,
               Fator = if_else(Condicao, "Dentro da circunferência", "Fora da circunferência"))

# Estimação do valor de pi
dados <- base %>% mutate(Acumulado = cumsum(Condicao),
                         piHat = 4*Acumulado/Order)

# Imprime-se valor estimado de pi na tela
cat("O valor estimado de pi foi:", round(dados$piHat[N], digits = 5))
## O valor estimado de pi foi: 3.116


Para uma melhor visualização do comportamento da estimação, pode-se gerar o gráfico de como se comporta o valor estimado de \(\pi\) a cada contabilização de pontos do vetor gerado (em azul). A linha vermelha representa o verdadeiro valor de \(\pi\):





Além disso, também pode ser gerado um diagrama para uma melhor visualização de como a simulação de Monte Carlo funciona:


Conclusão


Como esperado, os resultados encontrados foram bastante satisfatórios.

No exemplo exibido, o código com função de loop apresentou uma estimativa melhor que o código que utilizou vetores, mas isso se deve exclusivamente à semente do gerador aleatório. Em relação à performance computacional, o código vetorial apresenta tempo de execução bem mais reduzido.

Para uma melhor estimação, o valor do tamanho da amostra \(N\) pode ser aumentado, já que a estimação é gerada com base em uma convergência assintótica.

Mostra-se, assim, como é possível chegar em boas estimações de \(\pi\) através de métodos computacionais, em particular, o método de Monte Carlo.