간단한 몬테카를로 방법(Montecarlo Method) 예시

몬테카를로 방법은 과학과 공학 전 영역에 걸쳐 널리 사용되는 방법이다. 확률적인 해석을 요구하는 문제를 풀고싶을 때, 이 방법이 주로 쓰인다. 역사적으로는 물리학에서 자주 사용되었다고 전해진다. 1940년도에 원자로의 연쇄 반응 제어를 최초로 실현한 물리학자 엔리코 페르미는 중성자 확산(Neutron Diffusion)을 연구하면서, 이 방법을 자주 사용하였고, 로스앨러모스에서 맨하탄 프로젝트가 수행될 때 현대적인 버전의 몬테카를로 방법이 개발되었다고 전해진다.

참고 Montecarlo

몬테카를로 방법의 가장 간단한 예시로는 원주율(π)의 값을 추정하는 것이다.

넓이가 1인 정사각형을 생각하자. 정사각형의 한 꼭지점을 중심으로 반지름이 1인 사분원을 정사각형 안에 그린다. 그러면 사분원이 차지하는 넓이는 π/4가 될 것이다. 이제, 0 <= x <= 1인 x를 임의로 가져오고, 독립적으로 0 <= y <= 1인 y를 임의로 가져온 후, x^2 + y^2 <= 1일 확률은 사분원이 차지하는 넓이와 같은 값인 π/4가 된다.

이 과정을 여러 번 수행하는 알고리즘을 작성하고, 원주율의 값을 추정하시오.

Circle

※ 상대에게 상처를 주기보다 서로에게 도움이 될 수 있는 댓글을 달아 주세요.

35개의 풀이가 있습니다. 1 / 4 Page

import random
n=int(input("반복 횟수"))
count=0
for i in range(n):
    x=random.uniform(0,1.0)
    y=random.uniform(0,1.0)
    if (x**2+y**2)<=1:
        count=count+1
print("파이값은? ",4*count/n)
※ 상대에게 상처를 주기보다 서로에게 도움이 될 수 있는 댓글을 달아 주세요.
double pi=0;
@Test
public void 몬테카를로() {
    LongStream.rangeClosed(1l, 100l).forEach(v -> {
        pi=pi+BigDecimal.valueOf(LongStream.rangeClosed(0, v*3000).filter(x -> (Math.pow(Math.random(), 2) + Math.pow(Math.random(), 2)) <= 1).count()).divide(BigDecimal.valueOf(v*3000), MathContext.DECIMAL128).multiply(BigDecimal.valueOf(4)).doubleValue();
        logger.debug("{} : {}", v*3000, pi / v);
    });
}
@Test
public void 원주율추정() {
    double pi = 0;
    for (int i = 1; i < 100; i++) {
        pi = pi + reckoningPi(i * 3000);
        logger.debug("{} : {}", i * 3000, pi / i);
    }
}

private double reckoningPi(int limited) {
    double count = 0;
    for (double c = 0; c < limited; c++) {
        if (Math.pow(Math.random(), 2) + Math.pow(Math.random(), 2) <= 1) {
            count++;
        }
    }
    return count / limited * 4;
}
※ 상대에게 상처를 주기보다 서로에게 도움이 될 수 있는 댓글을 달아 주세요.

Ruby

m_pi = ->n { (1..n).reduce {|_,e|_+(rand(0..1.0)**2+rand(0..1.0)**2>1? 0:1) }*4.0/n }

Test

expect(m_pi[10000]).to be_within(0.1).of(3.1415)
expect(m_pi[20000]).to be_within(0.1).of(3.1415)
expect(m_pi[30000]).to be_within(0.1).of(3.1415)

Out

#=> p m_pi[30000]
3.1489333333333334
functional programming 에 대해 많이 배웁니다. 감사해요 - funnystyle, 2017/06/22 13:10 M D
※ 상대에게 상처를 주기보다 서로에게 도움이 될 수 있는 댓글을 달아 주세요.
#파이썬3.5.1
from random import *
from decimal import *

def montecarlo(n):
    points = []
    for i in range(n):
        points.append((Decimal(uniform(0,1)),Decimal(uniform(0,1))))
    c = 0
    for p in points:
        if p[0]**2 + p[1]**2 <= Decimal(1):
            c += 1
    return Decimal(c)/Decimal(n) * Decimal(4)

while True:
    n = input('How many repeat this program(If you want to stop, input nothing). : ')
    if n == '':
        break
    print('ㅠ will be',montecarlo(int(n))) # 파이를 ㅠ로 표기
※ 상대에게 상처를 주기보다 서로에게 도움이 될 수 있는 댓글을 달아 주세요.
public double getProbability(int repeat){
        int total = 0, count = 0;
        double x=0.0, y=0.0 , result = 0.0;
        total =  repeat;


        while(repeat > 0 ){
            repeat--;
            x = (double) Math.random()* 1.0;
            y = (double) Math.random()* 1.0;
            result = x * x + y * y;

            if(result <= 1.0){
                count++;
            }
            System.out.println("x : " + x +"," + "y:"+ y + " = " +result + "count : " + count);



        }
        return (double)count / (double) total;
    }

2016/05/14 15:07

xeo

※ 상대에게 상처를 주기보다 서로에게 도움이 될 수 있는 댓글을 달아 주세요.
import random
import math

N=100000
N_in_C=0

for i in range(N):
    x=random.random()
    y=random.random()
    if (x**2+y**2)<1 : N_in_C+=1

Pi_cal=((N_in_C/N)*4)
error=(abs(Pi_cal-math.pi)/math.pi)*100
print("Pi is %f" % Pi_cal)
print("Error rate = %f%%" % error)

※ 상대에게 상처를 주기보다 서로에게 도움이 될 수 있는 댓글을 달아 주세요.
dic = {}
for x in range(0, 1, 0.0001)
    for y in range(0, 1, 0.0001)
        if x**2 + y**2 =< 1:
            dic.append{x:y}
len(dic)/10**8 = pi/4
print(pi)

파이썬 처음인데 이게 맞을까요? ㅠㅠ

+1 DEPTH = 1000 i = 0 for x in range(0, DEPTH): for y in range(0, DEPTH): if x**2 + y**2 <= DEPTH**2: i += 1 pi = 4*i/(DEPTH**2) print(pi) 이걸 생각하신 듯 하네요. 근데 range() 함수는 소숫점 단위로 증가가 되지 않더라구요. 그래서 양변에 DEPTH의 제곱을 곱했다고 가정했습니다. 정밀도가 높아지면 확실히 괜찮겠지만, 문제에서는 확률적인 방법을 요구했네요. random 함수같은 것을 이용해 보세요. - Flair Sizz, 2016/06/17 23:17 M D
※ 상대에게 상처를 주기보다 서로에게 도움이 될 수 있는 댓글을 달아 주세요.
#include <iostream>
#include <stdlib.h>
#include <time.h>
#include <math.h>
using namespace std;
#define PI 3.141592
int main()
{
    srand((unsigned int)time(NULL));
    int n;
    float x, y;
    float mypi;
    float count = 0;

    cout << "몇 번?"; cin >> n;

    for (int i = 1; i <= n; i++)
    {
        x = (float)rand() / RAND_MAX;
        y = (float)rand() / RAND_MAX;
        if (pow(x, 2) + pow(y, 2) <= 1)
        {
            count++;
        }
    }
    mypi = count / n;
    cout << mypi * 4<< endl;

    return 0;
}
※ 상대에게 상처를 주기보다 서로에게 도움이 될 수 있는 댓글을 달아 주세요.
#include <iostream>
#include <stdlib.h>
#include <time.h>
#include <math.h>
using namespace std;
#define PI 3.141592
int main()
{
    srand((unsigned int)time(NULL));
    int n;
    float x, y;
    float mypi;
    float count = 0;

    cout << "몇 번?"; cin >> n;

    for (int i = 1; i <= n; i++)
    {
        x = (float)rand() / RAND_MAX;
        y = (float)rand() / RAND_MAX;
        if (pow(x, 2) + pow(y, 2) <= 1)
        {
            count++;
        }
    }
    mypi = count / n;
    cout << mypi * 4<< endl;

    return 0;
}
※ 상대에게 상처를 주기보다 서로에게 도움이 될 수 있는 댓글을 달아 주세요.
import random

count = 10000

def getPi(count):
    undercnt = 0
    for cnt in range(count):
        x = random.uniform(0,1.0)
        y = random.uniform(0,1.0)

        if x**2 + y**2 <= 1:
            undercnt +=1

    return (float(undercnt) / float(count))*4

if __name__ == '__main__':
    print getPi(count)

python 2.7

※ 상대에게 상처를 주기보다 서로에게 도움이 될 수 있는 댓글을 달아 주세요.

풀이 작성

※ 풀이작성 안내
  • 본문에 코드를 삽입할 경우 에디터 우측 상단의 "코드삽입" 버튼을 이용 해 주세요.
  • 마크다운 문법으로 본문을 작성 해 주세요.
  • 풀이를 읽는 사람들을 위하여 풀이에 대한 설명도 부탁드려요. (아이디어나 사용한 알고리즘 또는 참고한 자료등)
  • 작성한 풀이는 다른 사람(빨간띠 이상)에 의해서 내용이 개선될 수 있습니다.
목록으로
코딩도장

코딩도장은 프로그래밍 문제풀이를 통해서 코딩 실력을 수련(Practice)하는 곳입니다.


언어별 풀이 현황
전 체 x 35
python x 22
java x 3
ruby x 1
cpp x 5
go x 1
cs x 1
matlab x 1
javascript x 1