fork download
copy 
  1. #ifndef FIXEDPOINT_H
  2. #define FIXEDPOINT_H
  3.  
  4. //http://e...content-available-to-author-only...a.org/wiki/fixedb-point_arithmetic
  5.  
  6. //fixedb uses a specified number of bits for the integer part, and an underlying type.
  7. // if you want a specific number of bits (7) for the fractional part, use
  8. //  "BitsIn<unsigned>::Size-7" for the integer part
  9.  
  10. //fixed uses a specified number as the denominator for the fractional part
  11.  
  12. //fixedb is faster, but fixed gives a bigger range if the fractional part isn't 
  13. // multiples of powers of two, and rounds at more predictable places.
  14. //fixedb will round to the closest 128th of a dollar instead of the 100th.
  15.  
  16. //floating point only used as parameters to constructors (if desired) and in stream output
  17. //I intend to get accurate output later, eventually
  18.  
  19. #include <cassert>
  20. #include <iomanip>
  21. #include <iostream>
  22. #include <climits>
  23. #include <limits>
  24. #include <locale>
  25. #include <sstream>
  26. #include <type_traits>
  27.  
  28. #ifndef NFIXEDPOINTBOUNDS  //both NFIXEDPOINTBOUNDS and NDEBUG must be defined to get these asserts
  29. #define BOUNDASSERT(x)
  30. #else
  31. #define BOUNDASSERT(x) assert(x)
  32. #endif
  33.  
  34. //compile time logarithms
  35. template<unsigned long long num, unsigned long long base> struct ctlog { enum {value = ctlog<num/base, base>::value+1}; };
  36. template<unsigned long long base> struct ctlog<0, base> { enum {value = 0}; };
  37. //compile time exponents
  38. template<unsigned long long num, unsigned long long exp> struct ctpow { enum {value = num*ctpow<num-1, exp>::value}; };
  39. template<unsigned long long exp> struct ctpow<0, exp> { enum {value = 1}; };
  40.  
  41. //general case of maxvalue relies on signed rollover which is Undefined Behavior.  Please specialize
  42. template<class underlying, underlying one=1, bool biggerone=true> struct maxvalue { static const underlying value = maxvalue<underlying, one*2, (one*2>1)>::value; };
  43. template<class underlying, underlying one> struct maxvalue<underlying, one, false> { static const underlying value = one-underlying(1); };
  44. template<long long one>			struct maxvalue<long long, one, true> { static const long long value = LLONG_MAX; };
  45. template<unsigned long long one> struct maxvalue<unsigned long long, one, true> { static const unsigned long long value = ULLONG_MAX; };
  46. template<long one>				struct maxvalue<long, one, true> { static const long value = LONG_MAX; };
  47. template<unsigned long one>		struct maxvalue<unsigned long, one, true> { static const unsigned long value = ULONG_MAX; };
  48. template<int one>				struct maxvalue<int, one, true> { static const int value = INT_MAX; };
  49. template<unsigned one>			struct maxvalue<unsigned, one, true> { static const unsigned value = UINT_MAX; };
  50. template<short one>				struct maxvalue<short, one, true> { static const short value = SHRT_MAX; };
  51. template<unsigned short one>	struct maxvalue<unsigned short, one, true> { static const unsigned short value = USHRT_MAX; };
  52. template<char one>				struct maxvalue<char, one, true> { static const char value = CHAR_MAX; };
  53. template<signed char one>		struct maxvalue<signed char, one, true> { static const signed char value = SCHAR_MAX; };
  54. template<unsigned char one>		struct maxvalue<unsigned char, one, true> { static const unsigned char value = UCHAR_MAX; };
  55.  
  56. //used so I can pass values and integeral_const to the same function
  57. template<class underlying>
  58. struct carrier {
  59. 	underlying value;
  60. 	carrier(underlying v) :value(v) {}
  61. 	operator underlying() {return value;}
  62. };
  63.  
  64. // I KNOW I KNOW SHUT UP ALREADY COMPILER
  65. #pragma warning(push)
  66. #pragma warning(disable:4100) //warning C4100: 'times' : unreferenced formal parameter
  67. #pragma warning(disable:4127) //warning C4127: conditional expression is constant
  68. #pragma warning(disable:4146) //warning C4146: unary minus operator applied to unsigned type, result still unsigned
  69. #pragma warning(disable:4365) //warning C4365: 'initializing' : conversion from [STUFF] to 'const ulhstype', signed/unsigned mismatch
  70. #pragma warning(disable:4514) //warning C4514: [STUFF] : unreferenced inline function has been removed
  71. #pragma warning(disable:4710) //warning C4710: [STUFF] : function not inlined  (this pragma doesn't work)
  72. #pragma warning(disable:4723) //warning C4723: potential divide by 0
  73. #pragma warning(disable:4724) //warning C4724: potential mod by 0
  74. #pragma GCC diagnostic ignored "-Wtype-limits"
  75.  
  76. //multiplication
  77. //ignores signs, lhstype must be unsigned
  78. template<class lhstype, class timestype> 
  79. std::pair<lhstype, lhstype> longmultiply(lhstype lhs, timestype times) {
  80. 	static const lhstype bits = ctlog<maxvalue<lhstype>::value, 2>::value;
  81. 	static const lhstype half = bits/2+bits%2;
  82. 	static const lhstype mask = lhstype(1<<half)-1;
  83.  
  84. 	const lhstype umul((lhstype)(times.value>=0?times.value:-times.value)); 
  85. 	//               a      b      P1
  86.     //        *      c      d
  87. 	//          ----------------
  88. 	//             e(a*d) f(b*d)   P2
  89. 	//      g(a*c) h(b*c) 
  90. 	// -------------------------
  91. 	// i      j      k      l      P3
  92. 	//     m             n         P4
  93.  
  94. 	const lhstype a(lhs>>half); //P1 
  95. 	const lhstype b(lhs&mask);
  96. 	const lhstype c(umul>>half);
  97. 	const lhstype d(umul&mask);
  98. 	const lhstype e(a*d); //P2
  99. 	const lhstype f(b*d);
  100. 	const lhstype g(a*c);
  101. 	const lhstype h(b*c);
  102. 	const lhstype l(f&mask); //P3
  103. 	const lhstype k((e&mask) + (h&mask) + (f>>half));
  104. 	const lhstype j((g&mask) + (e>>half) + (h>>half) + (k>>half));
  105. 	const lhstype i((g>>half) + (j>>half));
  106. 	const lhstype m(((i&mask)<<half) | (j&mask)); //P4
  107. 	const lhstype n(((k&mask)<<half) | l);
  108. 	return std::pair<lhstype, lhstype>(m, n);
  109. }
  110. //mulitiplication specializations are to be preferred
  111. template<class timestype> 
  112. std::pair<unsigned int, unsigned int> longmultiply(unsigned int lhs, timestype times) {
  113. 	static_assert(sizeof(unsigned long long)>=sizeof(unsigned int)*2, "invalid assumption, comment out this overload");
  114. 	typedef unsigned long long unsignedlonglong; //Because GCC is stupid with function-style casts
  115. 	static const unsigned long long bits = CHAR_BIT*sizeof(unsigned int);
  116. 	unsigned long long r = unsignedlonglong(lhs) * unsignedlonglong(times.value);
  117. 	return std::pair<unsigned int, unsigned int>(r>>bits, r&((1ull<<bits)-1ull));
  118. }
  119.  
  120. template<class timestype> 
  121. std::pair<unsigned short, unsigned short> longmultiply(unsigned short lhs, timestype times) {
  122. 	static_assert(sizeof(unsigned int)>=sizeof(unsigned short)*2, "invalid assumption, comment out this overload");
  123. 	static const unsigned bits = CHAR_BIT*sizeof(unsigned short);
  124. 	unsigned int r = unsigned(lhs) * unsigned(times.value);
  125. 	return std::pair<unsigned short, unsigned short>(r>>bits, r&((1<<bits)-1));
  126. }
  127.  
  128. template<class timestype> 
  129. std::pair<unsigned char, unsigned char> longmultiply(unsigned char lhs, timestype times) {
  130. 	static_assert(sizeof(unsigned int)>=sizeof(unsigned char)*2, "invalid assumption, comment out this overload");
  131. 	static const unsigned bits = CHAR_BIT;
  132. 	unsigned int r = unsigned(lhs) * unsigned(times.value);
  133. 	return std::pair<unsigned char, unsigned char>(r>>bits, r&((1<<bits)-1));
  134. }
  135.  
  136. //division
  137. //ignores signs, lhstype must be unsigned
  138. template<class lhstype, class dividetype> 
  139. lhstype longdivision(std::pair<lhstype, lhstype> lhs, dividetype divide) {
  140. 	static const lhstype max = maxvalue<lhstype>::value;
  141. 	static const lhstype bits = ctlog<maxvalue<lhstype>::value, 2>::value;
  142. 	static const lhstype maxbits = max ^ (max>>2); // C000 0000
  143.  
  144. 	const lhstype udiv(divide.value>=0?divide.value:-divide.value); 
  145. 	BOUNDASSERT(lhs.first<udiv); //OVERFLOW
  146. 	lhstype result = 0;
  147. 	for(;;) {
  148. 		if (lhs.first == 0)
  149. 			return lhs.second/udiv+result;
  150. 		lhstype c = 0;
  151. 		int shift = 0;
  152. 		for(; ; shift+=2) {
  153. 			c = (lhs.first<<shift) | (shift?lhs.second>>(bits-shift):0);
  154. 			if (c&maxbits || shift==bits-2)
  155. 				break;
  156. 		}
  157. 		lhstype d = 0;
  158. 		if (udiv <= c) {
  159. 			d = (c/udiv);
  160. 			if (shift>0) // *= (max>>shift)+1 should result in *= 1, so...
  161. 				d *= (max>>shift)+1;
  162. 		} else if (shift<bits-2) {
  163. 			shift += 2;
  164. 			d = (c/(udiv/4));
  165. 			d *= (max>>shift)+1;
  166. 			assert(d>1);
  167. 			d -= 1; //sometimes it's one high?
  168. 		} else {
  169. 			d = 1;
  170. 		}
  171. 		assert(d>0);
  172. 		const std::pair<lhstype, lhstype> e(longmultiply(d, divide));
  173. 		assert((lhs.first > e.first) || (lhs.first==e.first && lhs.second>=e.second));
  174. 		result += d;
  175. 		if (e.second > lhs.second)
  176. 			lhs.first -= e.first+1;
  177. 		else
  178. 			lhs.first -= e.first;
  179. 		lhs.second -= e.second;
  180. 	}
  181.  }
  182.  
  183. //division specializations are to be preferred
  184. template<class dividetype> 
  185. unsigned int longdivision(std::pair<unsigned int, unsigned int> lhs, dividetype divide) {
  186. 	static_assert(sizeof(unsigned long long)>=sizeof(unsigned int)*2, "invalid assumption, comment out this overload");
  187. 	static const unsigned bits = CHAR_BIT*sizeof(unsigned int);
  188. 	typedef unsigned long long unsignedlonglong; //Because GCC is stupid with function-style casts
  189. 	unsigned long long t = unsignedlonglong(lhs.first) << bits | unsignedlonglong(lhs.second);
  190. 	return unsigned(t/divide.value);
  191. }
  192. template<class dividetype> 
  193. unsigned short longdivision(std::pair<unsigned short, unsigned short> lhs, dividetype divide) {
  194. 	static_assert(sizeof(unsigned int)>=sizeof(unsigned short)*2, "invalid assumption, comment out this overload");
  195. 	static const unsigned bits = CHAR_BIT*sizeof(unsigned short);
  196. 	unsigned int t = unsigned(lhs.first) << bits | unsigned(lhs.second);
  197. 	typedef unsigned short unsignedshort; //Because GCC is stupid with function-style casts
  198. 	return unsignedshort(t/divide.value);
  199. }
  200. template<class dividetype> 
  201. unsigned char longdivision(std::pair<unsigned char, unsigned char> lhs, dividetype divide) {
  202. 	static_assert(sizeof(unsigned int)>=sizeof(unsigned char)*2, "invalid assumption, comment out this overload");
  203. 	static const unsigned bits = CHAR_BIT*sizeof(unsigned char);
  204. 	unsigned int t = unsigned(lhs.first) << bits | unsigned(lhs.second);
  205. 	typedef unsigned char unsignedchar; //Because GCC is stupid with function-style casts
  206. 	return unsignedchar(t/divide.value);
  207. }
  208.  
  209. //do a multiplication and division
  210. template<class lhstype, class timestype, class dividetype>
  211. lhstype muldiv(lhstype lhs, timestype times, dividetype divide) {
  212. 	typedef typename std::make_unsigned<lhstype>::type ulhstype;
  213. 	const bool neg = (((lhs<0)!=(times.value<0))!=(divide.value<0));
  214. 	const std::pair<ulhstype, ulhstype> a = longmultiply<ulhstype, timestype>((ulhstype)(lhs>=0?lhs:-lhs), times);
  215. 	const lhstype b = (lhstype)longdivision(a, divide);
  216. 	return (neg ? -b : b);
  217. }	
  218.  
  219.     // base must be less than one tenth of the maximum value
  220. 	//base 16384 is: 
  221. 	// INT_MAX/16384 = 262143.999938965 maximum value
  222. 	// INT_MIN/16384 = 0 minimum value
  223. 	// lg(16384, 2) = 14 binary digits after point
  224. 	// 32-lg(16384, 2) = 18 binary digits above point
  225. 	// 14/3.322 = 4.21 decimal digits after point
  226. 	// 18/3.322 = 5.41 decimal digits above point (these should add to near 9.62)
  227. template<class base = std::integral_constant<unsigned, 16384> > 
  228. class fixedp
  229. {
  230. 	static_assert(base::value>0, "base must be greater than zero!");
  231. public:
  232. 	typedef typename base::value_type underlying;
  233. 	static const underlying fract = base::value;
  234. 	static const underlying max_value = maxvalue<underlying>::value;
  235. 	static const underlying min_value = std::numeric_limits<underlying>::is_signed ? -maxvalue<underlying>::value-1 : 0;
  236.     inline fixedp() :data(0) {}
  237. 	inline fixedp(const fixedp& rhs) :data(rhs.data) {}
  238.     inline fixedp(underlying rhs) : data(rhs*fract) {BOUNDASSERT(rhs.data<=max_value/fract);}
  239. #ifndef NO_HAS_DOUBLE
  240.     inline fixedp(double rhs) : data(underlying(rhs*fract+(1.0/fract/2.0))) {BOUNDASSERT(rhs<=max_value/fract && rhs>=min_value/fract);}
  241. #endif
  242. 	//second parameter ignored, used to signify that no math is required
  243.     inline fixedp(underlying numerator, underlying) : data(numerator) {} 
  244.     template<class rhsbase>
  245.     inline fixedp(const fixedp<rhsbase>& rhs) {assign<rhsbase, rhsbase::value/fract, rhsbase::value%fract, fract/rhsbase::value, fract%rhsbase::value>::asn(*this, rhs);}
  246.     inline ~fixedp() {}
  247.     inline fixedp& operator=(const fixedp& rhs) { data = rhs.data; return *this; }
  248.     inline fixedp  operator+(const fixedp& rhs) const { BOUNDASSERT((data>=0)?(max_value-data>=rhs.data):(rhs.data<min_value-data)); return fixedp(data+rhs.data,1); }
  249.     inline fixedp& operator+=(const fixedp& rhs) { BOUNDASSERT((data>=0)?(max_value-data>=rhs.data):(rhs.data<min_value-data)); data += rhs.data; return *this; }
  250.     inline fixedp  operator-(const fixedp& rhs) const { BOUNDASSERT(((data<0)==(rhs.data<0)) && (data>=0)?(max_value-data>=rhs.data):(min_value-data<rhs.data)); return fixedp(data-rhs.data,1); }
  251.     inline fixedp& operator-=(const fixedp& rhs) { BOUNDASSERT(((data<0)==(rhs.data<0)) && (data>=0)?(max_value-data>=rhs.data):(min_value-data<rhs.data)); data -= rhs.data; return *this; }
  252.     inline fixedp  operator*(underlying rhs) const { BOUNDASSERT(max_value/rhs.value>=data); return fixedp(data*rhs,1); }
  253.     inline fixedp& operator*=(underlying rhs)  { BOUNDASSERT(max_value/rhs.value>=data); data *= rhs; return *this; }
  254.     inline fixedp  operator*(const fixedp& rhs) const { return fixedp(muldiv(data, carrier<underlying>(rhs.data), base()), 1); }
  255. 	inline fixedp& operator*=(const fixedp& rhs) { data = muldiv(data, carrier<underlying>(rhs.data), base()); return *this;}
  256.     inline fixedp  operator/(underlying rhs) const {return fixedp(data/rhs,1); }
  257.     inline fixedp& operator/=(underlying rhs)  { data /= rhs; return *this; }
  258. 	inline fixedp  operator/(const fixedp& rhs) const {return fixedp(muldiv(data, base(), carrier<underlying>(rhs.data)), 1); }
  259. 	inline fixedp& operator/=(const fixedp& rhs) {data = muldiv(data, base(), carrier<underlying>(rhs.data)); return *this;}
  260.     inline fixedp  operator+() const { return fixedp(+data,1); }
  261.     inline fixedp  operator-() const { return fixedp(-data,1); }
  262.     inline fixedp& operator++() { BOUNDASSERT(max_value-fract>=data); data+=fract; return *this; }
  263.     inline fixedp  operator++(int) const { BOUNDASSERT(max_value-fract>=data); return fixedp(data+fract,1); }
  264.     inline fixedp& operator--() { BOUNDASSERT(min_value+fract<=data); data-=fract; return *this; }
  265.     inline fixedp  operator--(int) const { BOUNDASSERT(min_value+fract<=data); return fixedp(data-fract,1); }
  266.     inline bool    operator==(const fixedp& rhs) const { return data==rhs.data; }
  267.     inline bool    operator!=(const fixedp& rhs) const { return data!=rhs.data; }
  268.     inline bool    operator>(const fixedp& rhs) const { return data>rhs.data; }
  269.     inline bool    operator>=(const fixedp& rhs) const { return data>=rhs.data; }
  270.     inline bool    operator<=(const fixedp& rhs) const { return data<=rhs.data; }
  271.     inline bool    operator<(const fixedp& rhs) const { return data<rhs.data; }
  272.     inline bool    operator!() const { return !data; }
  273.     inline fixedp  operator~() const { return fixedp(~data,1); }
  274.     inline fixedp  operator&(const fixedp& rhs) const { return fixedp(data&rhs.data,1); }
  275.     inline fixedp& operator&=(const fixedp& rhs) { data &= rhs.data; return *this; }
  276.     inline fixedp  operator|(const fixedp& rhs) const { return fixedp(data|rhs.data,1); }
  277.     inline fixedp& operator|=(const fixedp& rhs) { data |= rhs.data; return *this; }
  278.     inline fixedp  operator^(const fixedp& rhs) const { return fixedp(data^rhs.data,1); }
  279.     inline fixedp& operator^=(const fixedp& rhs) { data ^= rhs.data; return *this; }
  280.     inline fixedp  operator<<(underlying rhs) const { BOUNDASSERT(rhs<CHAR_BIT*sizeof(underlying));BOUNDASSERT(max_value/(1<<rhs)>=data); return fixedp(data<<rhs,1); }
  281.     inline fixedp& operator<<=(underlying rhs) { BOUNDASSERT(rhs<CHAR_BIT*sizeof(underlying));BOUNDASSERT(max_value/(1<<rhs)>=data); data <<= rhs; return *this; }
  282.     inline fixedp  operator>>(underlying rhs) const { BOUNDASSERT(rhs<CHAR_BIT*sizeof(underlying)); return fixedp(data>>rhs,1); }
  283.     inline fixedp& operator>>=(underlying rhs) { BOUNDASSERT(rhs<CHAR_BIT*sizeof(underlying)); data >>= rhs; return *this; }
  284. 	//converts this type to any other type
  285.     template<typename convtype> convtype convert() const { return ((convtype)data)/((convtype)fract); }
  286. 	//gets the internal data
  287. 	inline underlying get_data() const {return data;}
  288. 	//not the best IO (actually it's pretty bad), but it works
  289.     template <class e, class t> friend std::basic_ostream<e,t>& operator<<(std::basic_ostream<e,t>& out, const fixedp& rhs) {
  290. 		out << rhs.data/fract;
  291. 		underlying left = (rhs.data%fract)*(rhs.data<0?-1:1);
  292.         if (left!=0 || out.flags()&std::ios_base::showpoint) {
  293.             out << '.';
  294. 			std::streamsize digits = out.precision();
  295. 			do {
  296. 				if (fract < max_value/2) {
  297. 					left*=10;
  298. 					out << char(left/fract+'0');
  299. 					left = left%fract;
  300. 				} else {
  301. 					underlying t = left/(fract/10);
  302. 					out << char(t+'0');
  303. 					left = left%(fract/10)*10;
  304. 				}
  305. 			} while(--digits>0 && left!=0);
  306. 		}
  307.         return out;
  308.     }
  309. 	//not the best IO (actually it's pretty bad), but it works
  310. 	template <class e, class t> friend std::basic_istream<e,t>& operator>>(std::basic_istream<e,t>& in, fixedp& rhs) {
  311. 		underlying s;
  312. 		in >> s;
  313. 		if (max_value/fract < s)
  314. 			in.setstate(std::basic_ios<e,t>::badbit | std::basic_ios<e,t>::failbit);
  315. 		else
  316. 			rhs.data = s*fract;
  317.         e next = e(in.peek());
  318. 		if (in && next=='.') {
  319. 			in.ignore(1);
  320. 			s=0;
  321. 			underlying pow10=1;
  322. 			for(next=e(in.peek()); in && std::isdigit(next, in.getloc()); next=e(in.peek())) {
  323. 				if (pow10>max_value/10) break;
  324. 				if (pow10>fract) break;
  325. 				s = s*(underlying)10+(underlying)(next-'0');
  326. 				pow10 *= (underlying)10;
  327. 				in.ignore(1);
  328. 			}
  329. 			rhs.data += muldiv(s, base(), carrier<underlying>(pow10)); //2677726520
  330. 			if (in && std::isdigit(next, in.getloc()) && next>='5')
  331. 				rhs.data += 1;
  332. 			for(; in && std::isdigit(next, in.getloc()); next=e(in.peek()))
  333. 				in.ignore(1);
  334. 		}
  335. 		return in;
  336.     }
  337. protected:
  338. 	//used for converting a fixedpoint of one base to another
  339. 	template<class rhsbase> friend class fixedp;
  340.     underlying data;
  341. 	template<class rhsbase, underlying div1, underlying mod1, underlying div2, underlying mod2>
  342. 	struct assign { static inline void asn(fixedp& left, fixedp<rhsbase> rhs) {left.data = muldiv(rhs.data, rhsbase::value, base::value);} };
  343. 	template<class rhsbase, underlying mod1, underlying div2>
  344. 	struct assign<rhsbase, 0, mod1, div2, 0> { static inline void asn(fixedp& left, fixedp<rhsbase> rhs)  {left.data = rhs.data/div2;} };
  345. 	template<class rhsbase, underlying div1, underlying mod2>
  346. 	struct assign<rhsbase, div1, 0, 0, mod2> { static inline void asn(fixedp& left, fixedp<rhsbase> rhs) {left.data = rhs.data*div1;} };
  347. 	template<class rhsbase>
  348. 	struct assign<rhsbase,1,0,1,0> { static inline void asn(fixedp& left, fixedp<rhsbase> rhs) {left.data = rhs.data;} };
  349. };
  350. #pragma warning(pop)
  351.  
  352. // // I would love to specialize this, but I heard I can't.  Commented out for now.
  353. //namespace std {
  354. //	template<class base>
  355. //	class numeric_limits<fixedp<base> > {
  356. //		typedef typename base::value_type underlying;
  357. //		typedef fixedp<base> myt;
  358. //	public:
  359. //		//static myt denorm_min( ) throw() {return myt(1,1);}
  360. //		static const int digits = numeric_limits<underlying>::digits;
  361. //		static const int digits10 = numeric_limits<underlying>::digits10;
  362. //		static myt epsilon() throw() {return myt(1,1);} 
  363. //		static const float_denorm_style has_denorm = denorm_absent;
  364. //		static const bool has_denorm_loss = false;
  365. //		static const bool has_infinity = false;
  366. //		static const bool has_quiet_NaN = false;
  367. //		static const bool has_signaling_NaN = false;
  368. //		//static myt infinity() throw() {return myt(numeric_limits<underlying>::infinity,1);}
  369. //		static const bool is_bounded = true;
  370. //		static const bool is_exact = true;
  371. //		static const bool is_iec559 = false;
  372. //		static const bool is_integer = false;
  373. //		static const bool is_modulo = true;
  374. //		static const bool is_signed = numeric_limits<underlying>::is_signed;
  375. //		static const bool is_specialized = true;
  376. //		static myt lowest() throw() {return myt(numeric_limits<underlying>::lowest(),1);}
  377. //		static myt max() throw() {return myt(numeric_limits<underlying>::max(),1);}
  378. //		static const int max_digits10 = numeric_limits<underlying>::max_digits10+2;
  379. //		static const int max_exponent = numeric_limits<underlying>::max_exponent/base::value;
  380. //		static const int max_exponent10 = numeric_limits<underlying>::max_exponent10/base::value;
  381. //		static myt min() throw() {return myt(numeric_limits<underlying>::min(),1);}
  382. //		static const int min_exponent = numeric_limits<underlying>::min_exponent/base::value;
  383. //		static const int min_exponent10 = numeric_limits<underlying>::min_exponent10/base::value;
  384. //		//static Type quiet_NaN( ) throw() {return max();}
  385. //		static const int radix = 2;
  386. //		static myt round_error() throw() {return myt(1,1);}
  387. //		static const float_round_style round_style = numeric_limits<underlying>::round_style;
  388. //		//static myt signaling_NaN() throw() {return max();}
  389. //		static const bool tinyness_before = false;
  390. //		static const bool traps = false;
  391. //	};
  392. //
  393. //	template<class base> struct is_arithmetic<fixedp<base>> : true_type {}
  394. //	template<class base> struct is_scalar<fixedp<base>> : true_type {}
  395. //	template<class base> struct is_signed<fixedp<base>> : is_signed<base::value_type>::value_type {}
  396. //	template<class base> struct is_unsigned<fixedp<base>> : is_unsigned<base::value_type>::value_type {}
  397. //	template<class <class, class> base, class underlying, underlying fract> 
  398. //	struct make_signed<fixedp<base<underlying, fract>>> {typename fixed<base<make_signed<base::value_type>::type, base::value> type;}
  399. //	template<class <class, class> base, class underlying, underlying fract> 
  400. //	struct make_unsigned<fixedp<base<underlying, fract>>> {typename fixed<base<make_unsigned<base::value_type>::type, base::value> type;}
  401. //}
  402.  
  403. #endif //FIXEDPOINT_H
  404.  
  405.  
  406.  
  407.  
  408.  
  409.  
  410.  
  411.  
  412.  
  413.  
  414.  
  415.  
  416.  
  417.  
  418.  
  419.  
  420.  
  421.  
  422.  
  423.  
  424.  
  425.  
  426.  
  427.  
  428.  
  429.  
  430.  
  431.  
  432.  
  433.  
  434.  
  435. #include <cassert>
  436. #include <iostream>
  437. #include <sstream>
  438. #include <ctime>
  439. #include <vector>
  440.  
  441. static const unsigned data_size = 50000;
  442. #ifdef _DEBUG
  443. #define NFIXEDPOINTBOUNDS 1
  444. static const unsigned max_test = 100000;
  445. #else
  446. static const unsigned max_test = 1000000;
  447. #endif
  448.  
  449. bool accuracy() {
  450. 	typedef fixedp<std::integral_constant<unsigned int, 2>> tinydenom;
  451. 	typedef fixedp<std::integral_constant<int, INT_MAX>> signedpercent; //2147483647
  452. 	typedef fixedp<std::integral_constant<unsigned int, UINT_MAX>> unsignedpercent; //4294967295	
  453. 	{
  454. 		std::cout << "unsignedpercent a = .1234567890" << std::endl;
  455. 		unsignedpercent a = .1234567890; //theory 0.1234567890
  456. 		assert(a.get_data() == 530242871u); //actual 0.12345678894780576229095458984375
  457. 		std::cout << "unsignedpercent b = 1.0" << std::endl;
  458. 		unsignedpercent b = 1.0; //theory 1
  459. 		assert(b.get_data() == 4294967295u); //actual 1.0
  460. 		std::cout << "unsignedpercent c = a*b" << std::endl;
  461. 		unsignedpercent c = a*b; //theory 0.1234567890
  462. 		assert(c.get_data() == 530242871u); //actual 0.12345678894780576229095458984375
  463. 		std::cout << "unsignedpercent d = .5" << std::endl;
  464. 		unsignedpercent d = .5; //theory 0.5
  465. 		assert(d.get_data() == 2147483647u); //actual 0.49999999988358467814596013122843
  466. 		std::cout << "unsignedpercent e = d*d" << std::endl;
  467. 		unsignedpercent e = d*d; //theory 0.25
  468. 		assert(e.get_data() == 1073741823u); //actual 0.24999999982537701721894019684264
  469. 		std::cout << "unsignedpercent f = a*a" << std::endl;
  470. 		unsignedpercent f = a*a; //theory 0.015241578750190521
  471. 		assert(f.get_data() == 65462082); //actual 0.015241578690531099841587967202437
  472. 	}
  473. 	{
  474. 		std::cout << "signedpercent a = .1234567890" << std::endl;
  475. 		signedpercent a = .1234567890; //theory .1234567890
  476. 		assert(a.get_data() == 265121435); //actual 0.12345678877246416582375027510512
  477. 		std::cout << "signedpercent b = 1.0" << std::endl;
  478. 		signedpercent b = 1.0; //theory 1.0
  479. 		assert(b.get_data() == 2147483647); //actual 1.0
  480. 		std::cout << "signedpercent c = a*b" << std::endl;
  481. 		signedpercent c = a*b; //theory 0.1234567890   (mul=07E6B74D 70329165, t=0FCD6E9B)
  482. 		assert(c.get_data() == 265121435); //actual 0.12345678877246416582375027510512
  483. 		std::cout << "signedpercent d = .5" << std::endl;
  484. 		signedpercent d = .5; //theory 0.5
  485. 		assert(d.get_data() == 1073741823); //actual 0.49999999976716935623771015379471
  486. 		std::cout << "signedpercent e = d*d" << std::endl;
  487. 		signedpercent e = d*d; //theory 0.25
  488. 		assert(e.get_data() == 536870911); //actual 0.24999999965075403435656523069207
  489. 		std::cout << "signedpercent f = -.5" << std::endl;
  490. 		signedpercent f = -.5; //theory -0.5
  491. 		assert(f.get_data() == -1073741823); //actual -0.49999999976716935623771015379471
  492. 		std::cout << "signedpercent g = e*f" << std::endl;
  493. 		signedpercent g = e*f; //theory -0.125
  494. 		assert(g.get_data() == -268435455); //actual -0.12499999959254637341599276914075
  495. 		std::cout << "signedpercent h = a*a" << std::endl;
  496. 		signedpercent h = a*a; //theory 0.015241578750190521
  497. 		assert(h.get_data() == 32731040); //actual 0.015241578228418518895478229455407
  498. 	}
  499. 	{
  500. 		std::cout << "tinydenom k = .5" << std::endl;
  501. 		tinydenom k = .5; //theory 0.5
  502. 		assert(k.get_data() == 1); //actual 0.5
  503. 		std::cout << "tinydenom l = k*k" << std::endl;
  504. 		tinydenom l = k*k; //theory 0.25
  505. 		assert(l.get_data() == 0); //actual 0.0
  506. 		std::cout << "tinydenom m = 2147483646u" << std::endl;
  507. 		tinydenom m = 2147483646u; //theory 2147483646.0
  508. 		assert(m.get_data() == 4294967292u); //actual 2147483646.0
  509. 		std::cout << "tinydenom n = k*m" << std::endl;
  510. 		tinydenom n = k*m; //theory 1073741823
  511. 		assert(n.get_data() == 2147483646u); //actual 1073741823.0
  512. 	}
  513. 	std::stringstream ss;
  514. 	std::string t;
  515. 	//float f(0);
  516. 	tinydenom td(0u);
  517. 	signedpercent sp(0);
  518. 	unsignedpercent up(0u);
  519.  
  520. 	ss.clear();
  521. 	ss << "2147483647.1234567890";
  522. 	std::cout << "ss >> td\n";
  523. 	ss >> td; 
  524. 	assert(td.get_data() == 4294967294u); //4294967295u => 2147483647
  525. 	ss.clear();
  526. 	ss << td;
  527. 	std::cout << "ss << td\n";
  528. 	ss >> t;
  529. 	assert(t=="2147483647");
  530.  
  531. 	ss.clear();
  532. 	ss << "0.1234567890";
  533. 	ss >> sp;
  534. 	assert(sp.get_data() == 265121435); //265121435 => 0.12345678877246416582375027510512
  535. 	ss.clear();
  536. 	ss << std::setprecision(9) << sp;
  537. 	ss >> t;
  538. 	assert(t=="0.123456789");
  539.  
  540. 	ss.clear();
  541. 	ss << "0.1234567890";
  542. 	std::cout << "ss >> up\n";
  543. 	ss >> up;
  544. 	assert(up.get_data() == 530242871u); // 530242871u => 0.12345678897655028593180474963314
  545. 	ss.clear();
  546. 	std::cout << "ss << up\n";
  547. 	ss << std::setprecision(9) << up;
  548. 	ss >> t;
  549. 	assert(t=="0.123456789");
  550.  
  551. 	return true;
  552. }
  553.  
  554. template<class Type, class integer>
  555. void Speed(const char* name) 
  556. {
  557. 	clock_t start, stop;
  558. 	Type total;
  559. 	std::cout << "timing " << name << std::endl;
  560.  
  561. 	srand(0);
  562. 	std::vector<Type> data(data_size);
  563. 	for(unsigned int i=0; i<data_size; ++i)
  564. 		data[i] += Type(double(rand())/128.0+1.0);
  565.  
  566. 	std::cout << "= ";
  567. 	total=integer(0);
  568. 	start = clock();
  569. 	for(unsigned i=0; i<max_test; ++i)
  570. 		total += data[i%data_size];
  571. 	stop = clock();
  572. 	std::cout << double(stop-start)/CLOCKS_PER_SEC << '\t' << total << std::endl;
  573.  
  574. 	std::cout << "+ ";
  575. 	total=integer(0);
  576. 	start = clock();
  577. 	for(unsigned i=0; i<max_test; ++i)
  578. 		total += data[i%data_size] + data[i*2%data_size];
  579. 	stop = clock();
  580. 	std::cout << double(stop-start)/CLOCKS_PER_SEC << '\t' << total << std::endl;
  581.  
  582. 	std::cout << "* ";
  583. 	total=integer(0);
  584. 	start = clock();
  585. 	for(unsigned i=0; i<max_test; ++i)
  586. 		total += data[i%data_size] * data[i*2%data_size];
  587. 	stop = clock();
  588. 	std::cout << double(stop-start)/CLOCKS_PER_SEC << '\t' << total << std::endl;
  589.  
  590. 	std::cout << "/ ";
  591. 	total=integer(0);
  592. 	start = clock();
  593. 	for(unsigned i=0; i<max_test; ++i)
  594. 		total += data[i%data_size] / data[i*2%data_size];
  595. 	stop = clock();
  596. 	std::cout << double(stop-start)/CLOCKS_PER_SEC << '\t' << total << std::endl;
  597. }
  598.  
  599.  
  600. int main() {
  601. 	typedef fixedp<std::integral_constant<char, 2>> smdenom;
  602. 	typedef fixedp<std::integral_constant<unsigned int, 11047>> prime;
  603. 	typedef fixedp<std::integral_constant<int, 65536>> pow2;
  604. 	typedef fixedp<std::integral_constant<int, 1000>> notpow2;
  605. 	typedef fixedp<std::integral_constant<unsigned int, 429496729>> bigdenom; //0x19999999
  606. 	typedef fixedp<std::integral_constant<short, 16>> slow2; //0x19999999
  607. 	typedef fixedp<std::integral_constant<short, 19>> slown2; //0x19999999
  608.  
  609. 	const bigdenom verify0(822083583, 1);                      //0822083583 0x30FFFFFF 1.9140625003456079871565215110171
  610. 	const bigdenom verify1(822083583, 1);                      //0822083583 0x30FFFFFF 1.9140625003456079871565215110171
  611. 	const bigdenom verify2 = verify0*verify1;                  //1573519358 0x5DC9FFFE 3.6636352543676764532472143693555
  612. 	const bigdenom verify3 = verify2/verify1;                  //0822083582 0x30FFFFFE 1.9140624980173015473652186999543
  613. 	std::cout << "0x822083583*0x822083583 = 3.663635254 = " << verify2 << std::endl;
  614. 	std::cout << "0x822083583/0x822083583 = 1.914062498 = " << verify3 << std::endl;
  615. 	std::cout << "UINT_MAX = " << UINT_MAX << " = " << ((unsigned)(UINT_MAX*1.0)) << std::endl;	
  616. 	std::cout << "INT_MAX = " << INT_MAX << " = " << maxvalue<int>::value << std::endl;
  617. 	std::cout << "signed bits = " << "31 = " << ctlog<maxvalue<int>::value, 2>::value << std::endl;
  618. 	std::cout << "half signed bits = " << "16 = " << ctlog<maxvalue<int>::value, 2>::value/2+ctlog<maxvalue<int>::value, 2>::value%2 << std::endl;
  619.  
  620. 	assert(accuracy());
  621.  
  622.     Speed<double, int>("double");
  623.     Speed<float, int>("float");
  624.     Speed<pow2, int>("int pow2");
  625.     Speed<notpow2, int>("int not2");
  626.     Speed<slow2, short>("short pow2");
  627.     Speed<slown2, short>("short not2");
  628.     return 0;
  629. }
  630.  
Success #stdin #stdout 0.35s 3180KB
comments (?)
stdin
copy 
Standard input is empty
stdout
copy 
0x822083583*0x822083583 = 3.663635254 = 3.663635
0x822083583/0x822083583 = 1.914062498 = 1.914062
UINT_MAX = 4294967295 = 4294967295
INT_MAX = 2147483647 = 2147483647
signed bits = 31 = 31
half signed bits = 16 = 16
unsignedpercent a = .1234567890
unsignedpercent b = 1.0
unsignedpercent c = a*b
unsignedpercent d = .5
unsignedpercent e = d*d
unsignedpercent f = a*a
signedpercent a = .1234567890
signedpercent b = 1.0
signedpercent c = a*b
signedpercent d = .5
signedpercent e = d*d
signedpercent f = -.5
signedpercent g = e*f
signedpercent h = a*a
tinydenom k = .5
tinydenom l = k*k
tinydenom m = 2147483646u
tinydenom n = k*m
ss >> td
ss << td
ss >> up
ss << up
timing double
= 0.01	8.3618e+12
+ 0	1.66974e+13
* 0.01	6.99028e+19
/ 0.01	6.72316e+06
timing float
= 0.01	8.3622e+12
+ 0	1.66975e+13
* 0.01	6.98868e+19
/ 0.01	6.72076e+06
timing int pow2
= 0	17113.625
+ 0.01	-30770.9375
* 0.02	30720
/ 0.02	-15050.056152
timing int not2
= 0	1571581.724
+ 0.01	128846.26
* 0.07	-1565070.788
/ 0.02	-1250827.54
timing short pow2
= 0.01	-1913.5
+ 0	1036.5
* 0.03	-1024
/ 0.03	-552.5
timing short not2
= 0	-863.578947
+ 0.01	514.105263
* 0.02	-1044.842105
/ 0.04	-393.684210
https://ideone.com/OkGql
language:
C++14 (gcc 8.3)
created:
12 years ago
visibility:
public

Share or Embed source code

Discover > Sphere Engine API

The brand new service which powers Ideone!

Discover > IDE Widget

Widget for compiling and running the source code in a web browser!