//
// Ldouble動作確認コード
//
#include <stdio.h>
#include <time.h>
#include <math.h>
#include "Ldouble.h"

void dump(char* str, double a, double b)
{
    printf("%s: %.16e: %.16e\n", str, a, (double)b);
}

int main()
{
    // 最初（又は必要な都度）Ldouble_initをCallすること
    Ldouble_init();

    double as = 10.;
    double bs = 123.;
    double cs;
    Ldouble al((unsigned __int64)10);
    Ldouble bl(bs);
    Ldouble cl;

    printf("as:= %.16e\nbs=: %.16e\n\n", as, bs);

    cs = as+bs; cl = al+bl;
    dump("+  ", cs, cl);
    cs = as-bs; cl = al-bl;
    dump("-  ", cs, cl);
    dump("abs", abs(cs), LD::abs(cl));
    if (as > bs) cs = bs - as;
    if (al > bl) cl = bl - al;
    dump("if ", cs, cl);
    cs = as*bs; cl = al*bl;
    dump("*  ", cs, cl);
    cs = bs/as; cl = bl/al;
    dump("/  ", cs, cl);
    cs = fmod(bs,as); cl = bl%al;
    dump("%  ", cs, cl);
    cs = 3.141592653589793238462643383279;
    cl = LD::pi();
    dump("pi ", cs, cl);
    dump("sin", sin(cs/6.0), LD::sin(cl/6.0));
    dump("cos", cos(cs/6.0), LD::cos(cl/6.0));
    LD::sincos(cl/6.0, al, bl);
    dump("sin", sin(cs/6.0), al);
    dump("cos", cos(cs/6.0), bl);
    dump("tan", tan(cs/6.0), LD::tan(cl/6.0));
    dump("ata", atan(cs/6.0), LD::atan(cl/6.0));
    dump("at2", atan2(cs/6.0, cs/9.0), LD::atan2(cl/6.0, cl/9.0));
    dump("sqr", sqrt(2.0), LD::sqrt(2.0));
    cs = exp(12.01);
    cl = LD::exp(12.01);
    dump("exp", cs, cl);
    dump("ln ", log(cs), LD::log(cl));
    cs = pow(cs, 0.1);
    cl = LD::pow(cl, 0.1);
    dump("pow", cs, cl);
    dump("int", int(cs), LD::rndint(cl));
    dump("dec", cs-int(cs), LD::rnddec(cl));

    printf ("\n");

    // マチンの公式(π)
    cs = 4.*(4.*atan(1./5.) - atan(1./239.));
    cl = 4.*(4.*LD::atan(Ldouble(1.)/5.) - LD::atan(Ldouble(1.)/239.));
    dump("Machin's pi", cs, cl);

    // マーダヴァ-ライプニッツのπ
    printf ("Leibniz's pi(5000000項での真値): 3.1415924535897932384646433832795\n");
    clock_t start,end;
    int i;
    start = clock();
    as=1.0; bs=2.0; cs=0.0;
    for (i=0; i<2500000; i++){
        cs += 1.0 / as;
        as += bs;
        cs -= 1.0 / as;
        as += bs;
    }
    cs *= 4;
    end = clock();
    printf ("Leibniz's pi(%d)  double(%.16e) time = %.3f[sec]\n", i*2, cs, (double)(end-start)/CLOCKS_PER_SEC);

    start = clock();
    al=1.0; bl=2.0; cl=0.0;
    for (i=0; i<2500000; i++){
        cl += 1.0 / al;
        al += bl;
        cl -= 1.0 / al;
        al += bl;
    }
    cl *= 4;
    end = clock();
    printf ("Leibniz's pi(%d) Ldouble(%.16e) time = %.3f[sec]\n", i*2, (double)cl, (double)(end-start)/CLOCKS_PER_SEC);

    return 0;
}