/***********************************************************************
Demonstrate some features of IEEE 754-2008 decimal floating-point
arithmetic.

Usage:
	dgcc decdemo.c -lmcw && ./a.out

[16-Aug-2010]
***********************************************************************/

#include <stdio.h>
#include <stdlib.h>
#include <mathcw.h>
#include <decfloat.h>

#define PRINTF	(void)printf

void
check_sizeof()
{
    if (sizeof(decimal_long_long_double) == sizeof(decimal_long_double))
	PRINTF("Warning: decimal_long_long_double and decimal_long_double "
	       "are identical %d-byte types on this platform\n\n", sizeof(decimal_long_double));
}

void
eps_df(void)
{
    decimal_float x;
    static const decimal_float ONE = 1.DF;
    static const decimal_float BASE = 1.e1DF;

    PRINTF("Demonstration of 32-bit machine-epsilon computation\n\n");

    x = ONE;

    while ((ONE + x / BASE) > ONE)
	x /= BASE;

    PRINTF("Computed decimal_float machine epsilon = %.7..3Hg\n", x);
    PRINTF("\n\n");
}

void
eps_dd(void)
{
    decimal_double x;
    static const decimal_double ONE = 1.DD;
    static const decimal_double BASE = 1.e1DD;

    PRINTF("Demonstration of 64-bit machine-epsilon computation\n\n");

    x = ONE;

    while ((ONE + x / BASE) > ONE)
	x /= BASE;

    PRINTF("Computed decimal_double machine epsilon = %.16..3DDg\n", x);
    PRINTF("\n\n");
}

void
eps_dl(void)
{
    decimal_long_double x;
    static const decimal_long_double ONE = 1.DL;
    static const decimal_long_double BASE = 1.e1DL;

    PRINTF("Demonstration of 128-bit machine-epsilon computation\n\n");

    x = ONE;

    while ((ONE + x / BASE) > ONE)
	x /= BASE;

    PRINTF("Computed decimal_long_double machine epsilon = %.34..3DLg\n", x);
    PRINTF("\n\n");
}

void
eps_dll(void)
{
    decimal_long_long_double x;
    static const decimal_long_long_double ONE = 1.DL;
    static const decimal_long_long_double BASE = 1.e1DL;

    PRINTF("Demonstration of 256-bit machine-epsilon computation\n\n");

    check_sizeof();

    x = ONE;

    while ((ONE + x / BASE) > ONE)
	x /= BASE;

    PRINTF("Computed decimal_long_long_double machine epsilon = %.70..3DLLg\n", x);
    PRINTF("\n\n");
}

void
limits_df(void)
{
    PRINTF("Decimal 32-bit floating-point limits\n\n");

    PRINTF("DEC_FLT_DEN        = %.7..3Hg\n", DEC_FLT_DEN);
    PRINTF("DEC_FLT_EPSILON    = %.7..3Hg\n", DEC_FLT_EPSILON);
    PRINTF("DEC_FLT_MANT_DIG   = %d\n",       DEC_FLT_MANT_DIG);
    PRINTF("DEC_FLT_MAX        = %.7..3Hg\n", DEC_FLT_MAX);
    PRINTF("DEC_FLT_MAX_EXP    = %d\n",       DEC_FLT_MAX_EXP);
    PRINTF("DEC_FLT_MIN        = %.7..3Hg\n", DEC_FLT_MIN);
    PRINTF("DEC_FLT_MIN_EXP    = %d\n",       DEC_FLT_MIN_EXP);
    PRINTF("\n\n");
}

void
limits_dd(void)
{
    PRINTF("Decimal 64-bit floating-point limits\n\n");

    PRINTF("DEC_DBL_DEN        = %.16..3DDg\n", DEC_DBL_DEN);
    PRINTF("DEC_DBL_EPSILON    = %.16..3DDg\n", DEC_DBL_EPSILON);
    PRINTF("DEC_DBL_MANT_DIG   = %d\n",         DEC_DBL_MANT_DIG);
    PRINTF("DEC_DBL_MAX        = %.16..3DDg\n", DEC_DBL_MAX);
    PRINTF("DEC_DBL_MAX_EXP    = %d\n",         DEC_DBL_MAX_EXP);
    PRINTF("DEC_DBL_MIN        = %.16..3DDg\n", DEC_DBL_MIN);
    PRINTF("DEC_DBL_MIN_EXP    = %d\n",         DEC_DBL_MIN_EXP);
    PRINTF("\n\n");
}

void
limits_dl(void)
{
    PRINTF("Decimal 128-bit floating-point limits\n\n");

    PRINTF("DEC_LDBL_DEN       = %.34..3DLg\n", DEC_LDBL_DEN);
    PRINTF("DEC_LDBL_EPSILON   = %.34..3DLg\n", DEC_LDBL_EPSILON);
    PRINTF("DEC_LDBL_MANT_DIG  = %d\n",         DEC_LDBL_MANT_DIG);
    PRINTF("DEC_LDBL_MAX       = %.34..3DLg\n", DEC_LDBL_MAX);
    PRINTF("DEC_LDBL_MAX_EXP   = %d\n",         DEC_LDBL_MAX_EXP);
    PRINTF("DEC_LDBL_MIN       = %.34..3DLg\n", DEC_LDBL_MIN);
    PRINTF("DEC_LDBL_MIN_EXP   = %d\n",         DEC_LDBL_MIN_EXP);
    PRINTF("\n\n");
}

void
limits_dll(void)
{
    PRINTF("Decimal 256-bit floating-point limits\n\n");

    check_sizeof();

    PRINTF("DEC_LLDBL_DEN      = %.70..3DLLg\n", DEC_LLDBL_DEN);
    PRINTF("DEC_LLDBL_EPSILON  = %.70..3DLLg\n", DEC_LLDBL_EPSILON);
    PRINTF("DEC_LLDBL_MANT_DIG = %d\n",          DEC_LLDBL_MANT_DIG);
    PRINTF("DEC_LLDBL_MAX      = %.70..3DLLg\n", DEC_LLDBL_MAX);
    PRINTF("DEC_LLDBL_MAX_EXP  = %d\n",          DEC_LLDBL_MAX_EXP);
    PRINTF("DEC_LLDBL_MIN      = %.70..3DLLg\n", DEC_LLDBL_MIN);
    PRINTF("DEC_LLDBL_MIN_EXP  = %d\n",          DEC_LLDBL_MIN_EXP);
    PRINTF("\n\n");
}

void
math_df(void)
{
    PRINTF("Demonstration of selected decimal 32-bit function values\n\n");

    PRINTF("cbrtdf  (22.DF)       = % .7..3Hg\n", cbrtdf(22.DF));
    PRINTF("elledf  (0.2DF)       = % .7..3Hg\n", elledf(0.2DF));
    PRINTF("erfcdf  (22.DF)       = % .7..3Hg\n", erfcdf(22.DF));
    PRINTF("pgammadf(10, 22.DF)   = % .7..3Hg\n", pgammadf(10, 22.DF));
    PRINTF("psidf   (22.DF)       = % .7..3Hg\n", psidf(22.DF));
    PRINTF("sbj0df  (22.DF)       = % .7..3Hg\n", sbj0df(22.DF));
    PRINTF("sindf   (22.DF)       = % .7..3Hg\n", sindf(22.DF));
    PRINTF("sqrtdf  (22.DF)       = % .7..3Hg\n", sqrtdf(22.DF));
    PRINTF("tgammadf(22.DF)       = % .7..3Hg\n", tgammadf(22.DF));
    PRINTF("\n\n");
}

void
math_dd(void)
{
    PRINTF("Demonstration of selected decimal 64-bit function values\n\n");

    PRINTF("cbrtd   (22.DD)       = % .16..3DDg\n", cbrtd(22.DD));
    PRINTF("elled   (0.2DD)       = % .16..3DDg\n", elled(0.2DD));
    PRINTF("erfcd   (22.DD)       = % .16..3DDg\n", erfcd(22.DD));
    PRINTF("pgammad (10, 22.DD)   = % .16..3DDg\n", pgammad(10, 22.DD));
    PRINTF("psid    (22.0F)       = % .16..3DDg\n", psid(22.DD));
    PRINTF("sbj0d   (22.DD)       = % .16..3DDg\n", sbj0d(22.DD));
    PRINTF("sind    (22.DD)       = % .16..3DDg\n", sind(22.DD));
    PRINTF("sqrtd   (22.DD)       = % .16..3DDg\n", sqrtd(22.DD));
    PRINTF("tgammad (22.DD)       = % .16..3DDg\n", tgammad(22.DD));
    PRINTF("\n\n");
}

void
math_dl(void)
{
    PRINTF("Demonstration of selected decimal 128-bit function values\n\n");

    PRINTF("cbrtdl  (22.DL)       = % .34..3DLg\n", cbrtdl(22.DL));
    PRINTF("elledl  (0.2DL)       = % .34..3DLg\n", elledl(0.2DL));
    PRINTF("erfcdl  (22.DL)       = % .34..3DLg\n", erfcdl(22.DL));
    PRINTF("pgammadl(10, 22.DL)   = % .34..3DLg\n", pgammadl(10, 22.DL));
    PRINTF("psidl   (22.0F)       = % .34..3DLg\n", psidl(22.DL));
    PRINTF("sbj0dl  (22.DL)       = % .34..3DLg\n", sbj0dl(22.DL));
    PRINTF("sindl   (22.DL)       = % .34..3DLg\n", sindl(22.DL));
    PRINTF("sqrtdl  (22.DL)       = % .34..3DLg\n", sqrtdl(22.DL));
    PRINTF("tgammadl(22.DL)       = % .34..3DLg\n", tgammadl(22.DL));
    PRINTF("\n\n");
}

void
math_dll(void)
{
    PRINTF("Demonstration of selected decimal 256-bit function values\n\n");

    check_sizeof();

    PRINTF("cbrtdll  (22.DLL)     = % .70..3DLLg\n", cbrtdll(22.DL));
    PRINTF("elledll  (0.2DLL)     = % .70..3DLLg\n", elledll(0.2DL));
    PRINTF("erfcdll  (22.DLL)     = % .70..3DLLg\n", erfcdll(22.DL));
    PRINTF("pgammadll(10, 22.DLL) = % .70..3DLLg\n", pgammadll(10, 22.DL));
    PRINTF("psidll   (22.0F)      = % .70..3DLLg\n",  psidll(22.DL));
    PRINTF("sbj0dll  (22.DLL)     = % .70..3DLLg\n", sbj0dll(22.DL));
    PRINTF("sindll   (22.DLL)     = % .70..3DLLg\n", sindll(22.DL));
    PRINTF("sqrtdll  (22.DLL)     = % .70..3DLLg\n", sqrtdll(22.DL));
    PRINTF("tgammadll(22.DLL)     = % .70..3DLLg\n", tgammadll(22.DL));
    PRINTF("\n\n");
}

void
math_hp(void)
{
    PRINTF("Expected values (from 200D symbolic-algebra evaluation)\n\n");

    PRINTF("cbrtdll  (22.DLL)     =  2.802_039_330_655_387_120_665_677_385_665_894_017_585_798_218_769_268_316_979_783_733_783_75\n");
    PRINTF("elledll  (0.2DLL)     =  1.554_968_546_242_529_283_474_427_007_626_648_437_452_982_275_568_916_084_827_105_580_383_232\n");
    PRINTF("erfcdll  (22.DLL)     =  1.621_905_860_933_472_513_052_034_647_026_123_265_177_675_596_121_851_615_118_740_214_899_267e-212\n");
    PRINTF("pgammadll(10, 22.DLL) = -1.702_523_439_631_811_782_184_094_486_897_304_774_841_852_577_408_045_311_669_762_286_963_616e-08\n");
    PRINTF("psidll   (22.0F)      =  3.068_143_039_861_196_669_924_876_026_445_032_981_842_169_957_058_100_049_372_679_705_171_979\n");
    PRINTF("sbj0dll  (22.DLL)     = -0.000_402_332_240_472_903_450_985_920_764_353_287_839_240_418_361_424_207_565_560_492_513_549_525_1\n");
    PRINTF("sindll   (22.DLL)     = -0.008_851_309_290_403_875_921_690_256_815_772_332_463_289_203_951_332_566_442_330_835_298_089_552\n");
    PRINTF("sqrtdll  (22.DLL)     =  4.690_415_759_823_429_554_565_630_113_544_466_280_588_228_353_411_737_153_605_701_891_017_025\n");
    PRINTF("tgammadll (22.DLL)    =  51_090_942_171_709_440_000\n");
    PRINTF("\n\n");
}

void
normalize_dl(void)
{
    PRINTF("Demonstration of normalization in decimal arithmetic\n\n");

    PRINTF("normalizedl(+0.00100DL)        = %s\n", ntosdl(normalizedl(+0.00100DL)));
    PRINTF("normalizedl(+1.00DL)           = %s\n", ntosdl(normalizedl(+1.00DL)));
    PRINTF("normalizedl(+100.DL)           = %s\n", ntosdl(normalizedl(+100.DL)));
    PRINTF("normalizedl(+100.00DL)         = %s\n", ntosdl(normalizedl(+100.00DL)));
    PRINTF("normalizedl(+qnandl(\"0x1234\")) = %s\n", ntosdl(normalizedl(+qnandl("0x1234"))));
    PRINTF("normalizedl(-qnandl(\"0x1234\")) = %s\n", ntosdl(normalizedl(-qnandl("0x1234"))));
    PRINTF("normalizedl(+inftydl())        = %s\n", ntosdl(normalizedl(+inftydl())));
    PRINTF("normalizedl(-inftydl())        = %s\n", ntosdl(normalizedl(-inftydl())));
    PRINTF("\n\n");
}

void
quantize_dl(void)
{
    PRINTF("Demonstration of quantization in decimal arithmetic\n\n");

    PRINTF("quantizedl(+1.DL, +1.00DL)                       = %s\n", ntosdl(quantizedl(+1.DL, +1.00DL)));
    PRINTF("quantizedl(+100.DL, +1.00DL)                     = %s\n", ntosdl(quantizedl(+100.DL, +1.00DL)));
    PRINTF("quantizedl(+0.125DL, +1.00DL)                    = %s\n", ntosdl(quantizedl(+0.125DL, +1.00DL)));
    PRINTF("quantizedl(+0.135DL, +1.00DL)                    = %s\n", ntosdl(quantizedl(+0.135DL, +1.00DL)));
    PRINTF("quantizedl(+0.145DL, +1.00DL)                    = %s\n", ntosdl(quantizedl(+0.145DL, +1.00DL)));
    PRINTF("quantizedl(+qnandl(\"0x1234\"), +1.00DL)           = %s\n", ntosdl(quantizedl(+qnandl("0x1234"), +1.00DL)));
    PRINTF("quantizedl(+100.DL, +qnandl(\"0x1234\"))           = %s\n", ntosdl(quantizedl(+100.DL, +qnandl("0x1234"))));
    PRINTF("quantizedl(+qnandl(\"0x1234\"), +qnandl(\"0x5678\")) = %s\n", ntosdl(quantizedl(+qnandl("0x1234"), +qnandl("0x5678"))));
    PRINTF("quantizedl(+inftydl(), -inftydl())               = %s\n", ntosdl(quantizedl(+inftydl(), -inftydl())));
    PRINTF("quantizedl(-inftydl(), +inftydl())               = %s\n", ntosdl(quantizedl(-inftydl(), +inftydl())));
    PRINTF("quantizedl(+100.00DL, -inftydl())                = %s\n", ntosdl(quantizedl(+100.00DL, -inftydl())));
    PRINTF("quantizedl(-inftydl(), +100.00DL)                = %s\n", ntosdl(quantizedl(-inftydl(), +100.00DL)));
    PRINTF("quantizedl(+1234567890123456.DL, +1.DL)          = %s\n", ntosdl(quantizedl(+1234567890123456.DL, +1.DL)));
    PRINTF("quantizedl(+1234567890123456.DL, +1.0DL)         = %s\n", ntosdl(quantizedl(+1234567890123456.DL, +1.0DL)));
    PRINTF("\n");

    PRINTF("samequantumdl(1.DL, 0.2DL)   = %d\n", samequantumdl(1.DL, 0.2DL));
    PRINTF("samequantumdl(1.DL, 2.DL)    = %d\n", samequantumdl(1.DL, 2.DL));
    PRINTF("samequantumdl(1.DL, 2.0DL)   = %d\n", samequantumdl(1.DL, 2.0DL));
    PRINTF("samequantumdl(1.DL, 2.00DL)  = %d\n", samequantumdl(1.DL, 2.00DL));
    PRINTF("samequantumdl(1.DL, 2.000DL) = %d\n", samequantumdl(1.DL, 2.000DL));
    PRINTF("\n");

    PRINTF("samequantumdl(+1.00DD, +9.99DD)                     = %d\n", samequantumdl(+1.00DD, +9.99DD));
    PRINTF("samequantumdl(+1.00DD, +9999.99DD)                  = %d\n", samequantumdl(+1.00DD, +9999.99DD));
    PRINTF("samequantumdl(+1.DD, +1.00DD)                       = %d\n", samequantumdl(+1.DD, +1.00DD));
    PRINTF("samequantumdl(+100.DD, +1.00DD)                     = %d\n", samequantumdl(+100.DD, +1.00DD));
    PRINTF("samequantumdl(+qnandl(\"0x1234\"), +1.00DD)           = %d\n", samequantumdl(+qnandl("0x1234"), +1.00DD));
    PRINTF("samequantumdl(+100.DD, +qnandl(\"0x1234\"))           = %d\n", samequantumdl(+100.DD, +qnandl("0x1234")));
    PRINTF("samequantumdl(+qnandl(\"0x1234\"), +qnandl(\"0x5678\")) = %d\n", samequantumdl(+qnandl("0x1234"), +qnandl("0x5678")));
    PRINTF("samequantumdl(+inftydl(), -inftydl())               = %d\n", samequantumdl(+inftydl(), -inftydl()));
    PRINTF("samequantumdl(-inftydl(), +inftydl())               = %d\n", samequantumdl(-inftydl(), +inftydl()));
    PRINTF("samequantumdl(+100.00DD, -inftydl())                = %d\n", samequantumdl(+100.00DD, -inftydl()));
    PRINTF("samequantumdl(-inftydl(), +100.00DD)                = %d\n", samequantumdl(-inftydl(), +100.00DD));
    PRINTF("\n\n");
}

void
sales_dl(void)
{
    static const decimal_long_double price_quantum = 0.01DL;
    static const volatile decimal_long_double item_cost = 0.70DL;
    static const volatile decimal_long_double sales_tax_rate = 0.05DL;
    volatile decimal_long_double sales_tax, total;


    sales_tax = item_cost * sales_tax_rate;;
    total = item_cost + sales_tax;

    PRINTF("Correct sales tax computation in decimal arithmetic\n\n");

    PRINTF("sales_tax_rate = %8.2DLf (printed)", sales_tax_rate);
    PRINTF(" = %s (stored)\n", ntosdl(sales_tax_rate));

    PRINTF("item_cost      = %8.2DLf (printed)", item_cost);
    PRINTF(" = %s (stored)", ntosdl(item_cost));
    PRINTF(" = %s (quantized)\n", ntosdl(quantizedl(item_cost, price_quantum)));

    PRINTF("sales_tax      = %8.2DLf (printed)", sales_tax);
    PRINTF(" = %s (computed)", ntosdl(sales_tax));
    PRINTF(" = %s (quantized)\n", ntosdl(quantizedl(sales_tax, price_quantum)));

    PRINTF("total          = %8.2DLf (printed)", total);
    PRINTF(" = %s (computed)", ntosdl(total));
    PRINTF(" = %s (quantized)\n", ntosdl(quantizedl(total, price_quantum)));
    PRINTF("\n\n");
}

int
main(void)
{
    limits_df();
    limits_dd();
    limits_dl();
    limits_dll();

    eps_df();
    eps_dd();
    eps_dl();
    eps_dll();

    math_df();
    math_dd();
    math_dl();
    math_dll();

    math_hp();

    sales_dl();

    normalize_dl();

    quantize_dl();

    return (EXIT_SUCCESS);
}