fnexp.hpp

浮点数基本运算 浮点数的基本运算主要有四则运算、符号处理、大小比较

// Exp and ln implementation -*- C++ -*-

// Copyright (C) 2005 Ben T. Bear
//
// This file is published under GNU GPL 2.0. See more in "copy.hpp".

// #include "copy.hpp"

#ifndef __bt2il_fn_exp_ln_hpp
#define __bt2il_fn_exp_ln_hpp

#include "fn.hpp"

namespace btil
{
  namespace fn
  {
    template <typename x>
    struct sqrt;

    template <typename x, int n>
    struct sqrt_n;

    template <typename x>
    struct exp;

    template <typename x>
    struct ln;

    template <typename x>
    struct lg;

    template <typename x, typename y>
    struct pow;

    template <typename x, int n>
    struct pow_n;

    template <int k, int n>
    struct pow_kn;



    template <int n>
    struct __factorial
    {
      static const double f_val = __factorial<n-1>::f_val * n;
    };

    template <>
    struct __factorial<0>
    {
      static const double f_val = 1.0;
    };



    template <typename a>
    struct __sqrt
    {
      static const int n2 = fn2em<a, 4>::e;
      static const double m = fn2em<a, 4>::m::f_val;

      static const double f0 = (m / 2.0 + 2.0) / 2;
      static const double f1 = (m / f0 + f0) / 2;
      static const double f2 = (m / f1 + f1) / 2;
      static const double f3 = (m / f2 + f2) / 2;
      static const double f4 = (m / f3 + f3) / 2;
      static const double f5 = (m / f4 + f4) / 2;
      static const double f6 = (m / f5 + f5) / 2;
      static const double f7 = (m / f6 + f6) / 2;
      static const double f = (m / f7 + f7) / 2;


      static const double f_val = f * pow_kn<2, n2>::value::f_val;
    };

    template <typename a>
    struct sqrt
    {
      struct value
      {
	static const double f_val = __sqrt<a>::f_val;
      };
    };



    template <int n>
    struct __e_item
    {
      static const double item = __e_item<n-1>::item / n;
      static const double e_n = __e_item<n-1>::e_n + item;
      //      static const double e_n = (__e_item<n-1>::e_n
      //			 + 1.0 / __factorial<n>::f_val);
    };

    template <>
    struct __e_item<0>
    {
      static const double item = 1.0;
      static const double e_n = __factorial<0>::f_val;
    };

    struct __e
    {
      static const double f_val = __e_item<22>::e_n;
    };



    template <typename x, int n>
    struct __exp_item
    {
      static const double item = __exp_item<x, n-1>::item * x::f_val / n;
      static const double exp_n = __exp_item<x, n-1>::exp_n + item;
    };

    template <typename x>
    struct __exp_item<x, 0>
    {
      static const double item = 1.0;
      static const double exp_n = item;
    };



    template <typename x, int n>
    struct __exp_sqr_p
    {
      static const double f_val = (__exp_sqr_p<x, n-1>::f_val
				   * __exp_sqr_p<x, n-1>::f_val);
    };
    template <typename x>
    struct __exp_sqr_p<x, 0>
    {
      static const double f_val = x::f_val;
    };
    template <typename x>
    struct sqrt;
    template <typename x, int n>
    struct __exp_sqr_m
    {
      struct __x
      {
	static const double f_val = __exp_sqr_m<x, n-1>::value::f_val;
      };
      struct value
      {
	static const double f_val = __sqrt<__x>::f_val;
      };
    };
    template <typename x>
    struct __exp_sqr_m<x, 0>
    {
      struct value
      {
	static const double f_val = x::f_val;
      };
    };
    template <bool neg, typename x, int n>
    struct __exp_sqr;
    template <typename x, int n>
    struct __exp_sqr<true, x, n>
    {
      static const double f_val = __exp_sqr_m<x, -n>::value::f_val;
    };
    template <typename x, int n>
    struct __exp_sqr<false, x, n>
    {
      static const double f_val = __exp_sqr_p<x, n>::f_val;
    };

    template <typename x>
    struct exp
    {
      struct __x
      {
	static const double __m = fn2em<x, 2>::m::f_val / 2;
	static const double f_val = x::f_val > 0 ? __m : -__m;
      };
      static const int __n = fn2em<x, 2>::e + 1;

      struct __a
      {
	static const double f_val = __exp_item<__x, 40>::exp_n;
      };

      struct value
      {
	static const double f_val = __exp_sqr<(__n<0), __a, __n>::f_val;
      };
    };


    
    struct __ln2;		// value of ln2

    template <int n>
    struct __ln2_item
    {
      static const double item = __ln2_item<n-1>::item / 9.0;
      static const double ln2_n = __ln2_item<n-1>::ln2_n + item / (2*n+1);
    };

    template <>
    struct __ln2_item<0>
    {
      static const double item = 2.0 / 3.0;
      static const double ln2_n = item;
    };

    struct __ln2
    {
      static const double f_val = __ln2_item<35>::ln2_n;
    };



    template <int n, typename x>
    struct __ln_item
    {
      static const double item = __ln_item<n-1, x>::item * -x::f_val;
      static const double ln_n = __ln_item<n-1, x>::ln_n + item / (n+1);
    };

    template <typename x>
    struct __ln_item<0, x>
    {
      static const double item = x::f_val;
      static const double ln_n = item;
    };

    template <typename x>
    struct ln
    {
      struct m
      {
	static const double f_val = fn2em<x, 2>::m::f_val / 2.0 - 1;
      };
      struct value
      {
	static const double f_val = (__ln_item<80, m>::ln_n
				     + __ln2::f_val * (fn2em<x, 2>::e+1));
      };
    };

    template <typename x, typename y>
    struct pow
    {
      struct __x
      {
	static const double f_val = y::f_val * ln<x>::value::f_val;
      };

      struct value
      {
	static const double f_val = exp<__x>::value::f_val;
      };
    };



    template <typename x, int n>
    struct __pow_n_p
    {
      static const double e = __pow_n_p<x, n/2>::f_val;
      static const double b = (n % 2 == 0) ? 1.0 : x::f_val;
      static const double f_val = e * e * b;
    };

    template <typename x>
    struct __pow_n_p<x, 0>
    {
      static const double f_val = 1.0;
    };

    template <bool neg, typename x, int n>
    struct __pow_n_bool;

    template <typename x, int n>
    struct __pow_n_bool<true, x, n>
    {
      static const double f_val = 1.0 / __pow_n_p<x, -n>::f_val;
    };

    template <typename x, int n>
    struct __pow_n_bool<false, x, n>
    {
      static const double f_val = __pow_n_p<x, n>::f_val;
    };

    template <typename x, int n>
    struct __pow_n
    {
      static const double f_val = __pow_n_bool<(n < 0), x, n>::f_val;
    };


    template <typename x, int n>
    struct pow_n
    {
      struct value
      {
	static const double f_val = __pow_n<x, n>::f_val;
      };
    };

    template <int k, int n>
    struct pow_kn
    {
      struct value
      {
	static const double f_val = __pow_kn<k, n>::value::f_val;
      };
    };



    template <typename x, int n>
    struct sqrt_n
    {
      struct __y
      {
	static const double f_val = 1.0 / n;
      };

      struct value
      {
	static const double f_val = pow<x, __y>::value::f_val;
      };
    };

    struct __ten
    {
      static const double f_val = 10.0;
    };
    struct __ln_ten
    {
      static const double f_val = ln<__ten>::value::f_val;
    };

    template <typename x>
    struct lg
    {
      struct value
      {
	static const double f_val = ln<x>::value::f_val / __ln_ten::f_val;
      };
    };
    struct __lg_e
    {
      static const double f_val = 1.0 / __ln_ten::f_val;
    };
  }
}

#endif