Métodos de Monte Carlo

Prévia do material em texto

Me´todos de Monte Carlo
Paulo Roberto de Carvalho Ju´nior
prcjunior@inf.ufpr.br
VRI – Visa˜o Robo´tica e Imagem
Universidade Federal do Parana´
Paulo Roberto de Carvalho Ju´nior Me´todos de Monte Carlo
Introduc¸a˜o
Me´todos de Monte Carlo sa˜o me´todos estat´ısticos baseados
em amostragens aleato´rias ou pseudoaleato´rias
Tempo de execuc¸a˜o determin´ıstico
Resultado aproximado
Quanto maior a quantidade de amostras, maior a
probabilidade de se aproximar do resultado correto
Paulo Roberto de Carvalho Ju´nior Me´todos de Monte Carlo
Exemplo: Ca´lculo de pi
Paulo Roberto de Carvalho Ju´nior Me´todos de Monte Carlo
Exemplo: Ca´lculo de pi
Equac¸a˜o da Circunfereˆncia: x2 + y2 = r2 ⇒ x2 + y2 = 1
A´rea do Quadrante do C´ırculo: AC =
pir2
4 =
pi
4
A´rea do Quadrado: AQ = r
2 = 1
AC
AQ
=
pi
4
1 =
pi
4 ⇒ pi = 4ACAQ
Paulo Roberto de Carvalho Ju´nior Me´todos de Monte Carlo
Algoritmo: Ca´lculo de pi
double calc_pi(int n)
{
int i;
int count = 0;
for (i = 0; i < n; ++i) {
double x = (double)rand() / RAND_MAX;
double y = (double)rand() / RAND_MAX;
if (x*x + y*y <= 1.0)
++ count;
}
return (double)count / n * 4.0;
}
Paulo Roberto de Carvalho Ju´nior Me´todos de Monte Carlo
Integrac¸a˜o pelo Me´todo de Monte Carlo
Problema: a complexidade da integrac¸a˜o nume´rica
determin´ıstica cresce exponencialmente de acordo com o
nu´mero de dimenso˜es do dom´ınio
Embora a integrac¸a˜o pelo me´todo de Monte Carlo tenha
resultado com erro probabil´ıstico, ela e´ escala´vel quanto a`
dimensa˜o do dom´ınio
Integrac¸a˜o pelo me´todo de Monte Carlo:∫
Ω
f (~x)d~x ≈ V 1
N
N∑
i=1
f (~xi )
onde Ω e´ um intervalo em Rm com volume V e ~xi sa˜o
amostras (pseudo)aleato´rias de Ω
Paulo Roberto de Carvalho Ju´nior Me´todos de Monte Carlo
Integrac¸a˜o pelo Me´todo de Monte Carlo
Paulo Roberto de Carvalho Ju´nior Me´todos de Monte Carlo
Algoritmo: Integrac¸a˜o pelo Me´todo de Monte Carlo
double calc_integral(double (*f)(double),
double a, double b, int n)
{
int i;
double fsum = 0.0;
for (i = 0; i < n; ++i) {
double x = a + (( double)rand()/RAND_MAX)
* (b - a);
fsum += f(x);
}
return (fsum / n) * (b - a);
}
Paulo Roberto de Carvalho Ju´nior Me´todos de Monte Carlo
Geradores de Nu´meros Pseudoaleato´rios
Sequeˆncias de nu´meros verdadeiramente aleato´rios:
Imprevis´ıveis
Irreprodut´ıveis
Sa˜o dif´ıceis de se obter
Requer hardware especializado
Taxa de transfereˆncia da sequeˆncia e´ geralmente baixa
Sequeˆncias de nu´meros pseudoaleato´rios:
Sa˜o geradas por algoritmos
Sa˜o determin´ısticas e reprodut´ıveis
Um mesmo valor de semente gera sempre a mesma sequeˆncia
Na˜o sa˜o aleato´rias, mas possuem propriedades pro´ximas de
sequeˆncias aleato´rias
Paulo Roberto de Carvalho Ju´nior Me´todos de Monte Carlo
Qualidade de Geradores de Nu´meros Pseudoaleato´rios
A qualidade de uma sequeˆncia pseudoaleato´ria e´ determinante
para a qualidade do resultado de um algoritmo de Monte Carlo
O per´ıodo da sequeˆncia deve ser suficientemente longo
Deve se comportar estatisticamente como a uma sequeˆncia de
nu´meros verdadeiramente aleato´rios
Deve estar equilibrada dentro de uma determinada distribuic¸a˜o
probabil´ıstica
Uniforme
Gaussiana (Normal)
Poisson
A relac¸a˜o entre os nu´meros da sequeˆncia deve ser dif´ıcil de se
identificar
Deve passar em testes emp´ıricos de aleatoriedade
O ca´lculo do pro´ximo nu´mero da sequeˆncia deve ser eficiente
e ra´pido
Paulo Roberto de Carvalho Ju´nior Me´todos de Monte Carlo
Gerador Congruencial Linear
Gerador de nu´meros pseudoaleato´rios de distribuic¸a˜o uniforme
Fa´cil de implementar
Eficiente e ra´pido
xn+1 = (a · xn + c) mod m
onde:
m, m > 0: mo´dulo
a, 0 < a < m: multiplicador
c , 0 6 c < m: incremento
x0, 0 6 x0 < m: semente
Paulo Roberto de Carvalho Ju´nior Me´todos de Monte Carlo
Gerador Congruencial Linear
A qualidade da sequeˆncia depende dos valores escolhidos para
as constantes
O per´ıodo da sequeˆncia e´ limitado por m, mas pode ser muito
menor se o valor do multiplicador a e do incremento c na˜o
forem escolhido adequadamente
De acordo com o teorema de Hull-Dobell, os seguintes
crite´rios devem ser obedecidos para maximizar o per´ıodo:
c e m devem ser primos entre si
a− 1 deve ser divis´ıvel por todos os fatores primos de m
a− 1 deve ser mu´ltiplo de 4 se m e´ mu´ltiplo de 4
Exemplos:
xn+1 = (3 · xn + 5) mod 32, x0 = 9
glibc: xn+1 = (1103515245 · xn + 12345) mod 231
Paulo Roberto de Carvalho Ju´nior Me´todos de Monte Carlo
Bibliotecas de Geradores de Nu´meros Pseudoaleato´rios
glibc (C/C++):
Per´ıodo: 231
#include <stdlib.h>
srand(): Inicializa a semente
rand(): Retorna um nu´mero pseudoaleato´rio entre 0 e
RAND MAX
SPRNG (C/C++, FORTRAN):
Inclui va´rios geradores diferentes, com diferentes per´ıodos
(2219, 264, 261, 248 e 231)
Suporte a processamento paralelo e concorrente
Dispon´ıvel em http://www.sprng.org/
Paulo Roberto de Carvalho Ju´nior Me´todos de Monte Carlo
Rotinas U´teis
Usar sementes diferente para cada vez que executar o
programa (resoluc¸a˜o de 1s)
srand(time(NULL));
var ← random[0, 1]
var = (double)rand() / RAND MAX;
var ← random[a, b]
var = a + (b - a) * ((double)rand() / RAND MAX);
Paulo Roberto de Carvalho Ju´nior Me´todos de Monte Carlo
Exerc´ıcio 1
Integre pelo me´todo de Monte Carlo a func¸a˜o f (x , y) definida no
intervalo [(−pi,−pi), (pi, pi)]:
f (x , y) = 3x2 + 5y2 + 3sin(x3) + 5cos(y3)
O programa deve utilizar 100 milho˜es de amostras de pontos
pseudoaleato´rios em f (x , y)
Use como semente para srand() o nu´mero 1234
Calcule o wall-clock time da integrac¸a˜o
Acesse http://www.inf.ufpr.br/prcjunior/ para baixar
o co´digo base
Paulo Roberto de Carvalho Ju´nior Me´todos de Monte Carlo
Exerc´ıcio 2
Otimize a func¸a˜o de integrac¸a˜o implementada no exerc´ıcio anterior
O co´digo pode ser otimizado especificamente para a func¸a˜o
f (x , y) apesentada, ou seja, na˜o e´ necessa´rio deixa´-lo gene´rico
para integrar qualquer func¸a˜o
Utilize as estrate´gias de otimizac¸a˜o aprendidas na disciplina
Utilize a ferramenta likwid-perfctr para auxiliar na ana´lise de
desempenho
Quem consegue o menor tempo de execuc¸a˜o? :)
Paulo Roberto de Carvalho Ju´nior Me´todos de Monte Carlo