테일러 급수를 이용한 sin() 및 BBP를 이용한 원주율 구현

이전 블로그에서 적었듯, ln()을 구현한 김에 sin()도 함께 구현해봤고, 원주율 계산도 간단하게 구현해봤다.

 

1. sin(x)의 구현

 

물론, 이번에도 테일러 급수다.

sin(x)는 아래와 같이 전개된다.

 

 

이번에도 역시 구현 자체는 그리 어렵지 않지만 적절한 횟수를 정하는 것이 필요하다.

 

const static double PI = 3.141592653589796323846;

void SIN_v1(double x, double limit, double& result, int& count) {
    result = 0;
    count = 0;

    while (x > 2 * PI) { x -= 2 * PI; }
    while (x < 0) { x += 2 * PI; }

    double xxxx = x * x * x * x;
    
    double t1 = x;
    double t2 = -x * x * x / 6;
    
    result = t1 + t2;
    count = 0;
    int i4 = 0;
    while (true) {
        ++count;
        i4 += 4;

        t1 = t1 * xxxx / ((i4 - 2) * (i4 - 1) * i4 * (i4 + 1));
        t2 = t2 * xxxx / (i4 * (i4 + 1) * (i4 + 2) * (i4 + 3));

        if (isinf(t1) || isinf(t2)) {
            break;
        }

        result += (t1 + t2);

        if (fabs(t1) < limit) {
            break;
        }
    }
}

 

ln()과 유사하게, 약 10회 정도 루프를 돌리니 상당히 정확한 값이 나온다.

이 점을 고려해서 작성한 결과는 아래와 같다.

 

const static double PI = 3.141592653589796323846;

double SIN_v2(double x) {
    while (x > 2 * PI) { x -= 2 * PI; }
    while (x < 0) { x += 2 * PI; }

    double xxxx = x * x * x * x;

    double t1 = x;
    double t2 = -x * x * x / 6;

    double result = t1 + t2;
    int i4 = 0;
    for (int i = 0; i < 20; ++i) {
        i4 += 4;

        t1 = t1 * xxxx / ((i4 - 2) * (i4 - 1) * i4 * (i4 + 1));
        t2 = t2 * xxxx / (i4 * (i4 + 1) * (i4 + 2) * (i4 + 3));

        if (isinf(t1) || isinf(t2)) {
            break;
        }

        result += (t1 + t2);
    }

    return result;
}

 

2. pi()의 구현

 

뉴튼과 같은 시대에 미적분학을 창시한 라이프니츠^[각주:1]는 원주율을 아래와 같이 전개하였다.

 

 

약 두 세대 쯤 뒤의 수학자인 그 위대하신 오일러^[각주:2]는 아래와 같이 전개하였다.

 

 

1995년에는 데이빗 베일리 - 피터 보어와인 - 시몽 플루프가 π에 관련된 새로운 무한급수를 발견^[각주:3]했다.

 

 

데이빗 베일리는 수학자이기도 했지만, 또한 나사에서 14년을 컴퓨터 과학자(Computer Scientist)로서도 근무했었다.

이런 그의 특기가 반영된 식답게 수렴속도도 빠르다.

 

루프의 적절한 횟수를 정하기 위해 일단 작성한 코드는 아래와 같다.

 

void PI_v1(double limit, double& result, int& count) {
    count = 0;

    result = 4 - 0.5 - 0.2 - (1 / 6.0);
    int sixteen = 1;

    double temp;
    int i8 = 0;
    while (true) {
        ++count;
        i8 += 8;
        sixteen *= 16;
        temp = ((4.0 / (i8 + 1)) - (2.0 / (i8 + 4)) - (1.0 / (i8 + 5)) - (1.0 / (i8 + 6))) / sixteen;

        if (isinf(temp)) {
            --count;
            break;
        }

        result += temp;

        if (temp < limit) {
            break;
        }
    }
}

 

이렇게 돌려보면 7회까지만 돌릴 수 있다.

이 점을 고려해서 작성한 코드는 아래와 같다.

 

double PI_v2() {
    double result = 4 - 0.5 - 0.2 - (1 / 6.0);
    int sixteen = 1;

    double temp;
    int i8 = 0;
    for (int i = 0; i < 7; ++i) {
        i8 += 8;
        sixteen *= 16;
        temp = ((4.0 / (i8 + 1)) - (2.0 / (i8 + 4)) - (1.0 / (i8 + 5)) - (1.0 / (i8 + 6))) / sixteen;

        result += temp;
    }

    return result;
}

 

이렇게 계산된 결과는 아래와 같다.

 

 

 

1. 미적분학 외에도 수많은 수학, 철학 분야에 지대한 업적이 있는 그야말로 시대를 뛰어넘는 핵천재. 미적분학 기호 dy/dx가 바로 라이프니츠 표기법 [본문으로]
2. 백내장으로 두 눈 모두를 실명한 상태에서도 실명 이전보다 더 많은 논문을 써낸 머신. 자신이 쓰고 읽은 모든 수학책을 다 외운, 더 이상의 설명이 불가능한 천재 오브 천재. 우리가 사용하는 대부분의 수학 기호는 오일러 표기법이라 봐도 무방함. 자연대수 e는 Euler의 앞글자이며, π 라는 표현도 오일러가 명명함 [본문으로]
3. 엄밀히는 이 무한급수는 단순히 원주율 뿐만 아니라 다양한 초월함수를 기술하는 급수임 [본문으로]