yum-mirror/VkFFT

git clone https://git.yummers.dev/yum-mirror/VkFFT

Tolmachev DmitriiBugfixes -Issue #50 should be partially resolved -Issue #51 should be resolved -Fixed ROCm 4.5 build problems with new version of HIPRTC -Updated Half version -Fixed some missing error handling and tempStr deallocation2ce8588

master
206.6 KiB4601 linesraw
1// half - IEEE 754-based half-precision floating-point library.
2//
3// Copyright (c) 2012-2021 Christian Rau <rauy@users.sourceforge.net>
4//
5// Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation 
6// files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, 
7// modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the 
8// Software is furnished to do so, subject to the following conditions:
9//
10// The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
11//
12// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE 
13// WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR 
14// COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, 
15// ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
16
17// Version 2.2.0
18
19/// \file
20/// Main header file for half-precision functionality.
21
22#ifndef HALF_HALF_HPP
23#define HALF_HALF_HPP
24
25#define HALF_GCC_VERSION (__GNUC__*100+__GNUC_MINOR__)
26
27#if defined(__INTEL_COMPILER)
28	#define HALF_ICC_VERSION __INTEL_COMPILER
29#elif defined(__ICC)
30	#define HALF_ICC_VERSION __ICC
31#elif defined(__ICL)
32	#define HALF_ICC_VERSION __ICL
33#else
34	#define HALF_ICC_VERSION 0
35#endif
36
37// check C++11 language features
38#if defined(__clang__)										// clang
39	#if __has_feature(cxx_static_assert) && !defined(HALF_ENABLE_CPP11_STATIC_ASSERT)
40		#define HALF_ENABLE_CPP11_STATIC_ASSERT 1
41	#endif
42	#if __has_feature(cxx_constexpr) && !defined(HALF_ENABLE_CPP11_CONSTEXPR)
43		#define HALF_ENABLE_CPP11_CONSTEXPR 1
44	#endif
45	#if __has_feature(cxx_noexcept) && !defined(HALF_ENABLE_CPP11_NOEXCEPT)
46		#define HALF_ENABLE_CPP11_NOEXCEPT 1
47	#endif
48	#if __has_feature(cxx_user_literals) && !defined(HALF_ENABLE_CPP11_USER_LITERALS)
49		#define HALF_ENABLE_CPP11_USER_LITERALS 1
50	#endif
51	#if __has_feature(cxx_thread_local) && !defined(HALF_ENABLE_CPP11_THREAD_LOCAL)
52		#define HALF_ENABLE_CPP11_THREAD_LOCAL 1
53	#endif
54	#if (defined(__GXX_EXPERIMENTAL_CXX0X__) || __cplusplus >= 201103L) && !defined(HALF_ENABLE_CPP11_LONG_LONG)
55		#define HALF_ENABLE_CPP11_LONG_LONG 1
56	#endif
57#elif HALF_ICC_VERSION && defined(__INTEL_CXX11_MODE__)		// Intel C++
58	#if HALF_ICC_VERSION >= 1500 && !defined(HALF_ENABLE_CPP11_THREAD_LOCAL)
59		#define HALF_ENABLE_CPP11_THREAD_LOCAL 1
60	#endif
61	#if HALF_ICC_VERSION >= 1500 && !defined(HALF_ENABLE_CPP11_USER_LITERALS)
62		#define HALF_ENABLE_CPP11_USER_LITERALS 1
63	#endif
64	#if HALF_ICC_VERSION >= 1400 && !defined(HALF_ENABLE_CPP11_CONSTEXPR)
65		#define HALF_ENABLE_CPP11_CONSTEXPR 1
66	#endif
67	#if HALF_ICC_VERSION >= 1400 && !defined(HALF_ENABLE_CPP11_NOEXCEPT)
68		#define HALF_ENABLE_CPP11_NOEXCEPT 1
69	#endif
70	#if HALF_ICC_VERSION >= 1110 && !defined(HALF_ENABLE_CPP11_STATIC_ASSERT)
71		#define HALF_ENABLE_CPP11_STATIC_ASSERT 1
72	#endif
73	#if HALF_ICC_VERSION >= 1110 && !defined(HALF_ENABLE_CPP11_LONG_LONG)
74		#define HALF_ENABLE_CPP11_LONG_LONG 1
75	#endif
76#elif defined(__GNUC__)										// gcc
77	#if defined(__GXX_EXPERIMENTAL_CXX0X__) || __cplusplus >= 201103L
78		#if HALF_GCC_VERSION >= 408 && !defined(HALF_ENABLE_CPP11_THREAD_LOCAL)
79			#define HALF_ENABLE_CPP11_THREAD_LOCAL 1
80		#endif
81		#if HALF_GCC_VERSION >= 407 && !defined(HALF_ENABLE_CPP11_USER_LITERALS)
82			#define HALF_ENABLE_CPP11_USER_LITERALS 1
83		#endif
84		#if HALF_GCC_VERSION >= 406 && !defined(HALF_ENABLE_CPP11_CONSTEXPR)
85			#define HALF_ENABLE_CPP11_CONSTEXPR 1
86		#endif
87		#if HALF_GCC_VERSION >= 406 && !defined(HALF_ENABLE_CPP11_NOEXCEPT)
88			#define HALF_ENABLE_CPP11_NOEXCEPT 1
89		#endif
90		#if HALF_GCC_VERSION >= 403 && !defined(HALF_ENABLE_CPP11_STATIC_ASSERT)
91			#define HALF_ENABLE_CPP11_STATIC_ASSERT 1
92		#endif
93		#if !defined(HALF_ENABLE_CPP11_LONG_LONG)
94			#define HALF_ENABLE_CPP11_LONG_LONG 1
95		#endif
96	#endif
97	#define HALF_TWOS_COMPLEMENT_INT 1
98#elif defined(_MSC_VER)										// Visual C++
99	#if _MSC_VER >= 1900 && !defined(HALF_ENABLE_CPP11_THREAD_LOCAL)
100		#define HALF_ENABLE_CPP11_THREAD_LOCAL 1
101	#endif
102	#if _MSC_VER >= 1900 && !defined(HALF_ENABLE_CPP11_USER_LITERALS)
103		#define HALF_ENABLE_CPP11_USER_LITERALS 1
104	#endif
105	#if _MSC_VER >= 1900 && !defined(HALF_ENABLE_CPP11_CONSTEXPR)
106		#define HALF_ENABLE_CPP11_CONSTEXPR 1
107	#endif
108	#if _MSC_VER >= 1900 && !defined(HALF_ENABLE_CPP11_NOEXCEPT)
109		#define HALF_ENABLE_CPP11_NOEXCEPT 1
110	#endif
111	#if _MSC_VER >= 1600 && !defined(HALF_ENABLE_CPP11_STATIC_ASSERT)
112		#define HALF_ENABLE_CPP11_STATIC_ASSERT 1
113	#endif
114	#if _MSC_VER >= 1310 && !defined(HALF_ENABLE_CPP11_LONG_LONG)
115		#define HALF_ENABLE_CPP11_LONG_LONG 1
116	#endif
117	#define HALF_TWOS_COMPLEMENT_INT 1
118	#define HALF_POP_WARNINGS 1
119	#pragma warning(push)
120	#pragma warning(disable : 4099 4127 4146)	//struct vs class, constant in if, negative unsigned
121#endif
122
123// check C++11 library features
124#include <utility>
125#if defined(_LIBCPP_VERSION)								// libc++
126	#if defined(__GXX_EXPERIMENTAL_CXX0X__) || __cplusplus >= 201103
127		#ifndef HALF_ENABLE_CPP11_TYPE_TRAITS
128			#define HALF_ENABLE_CPP11_TYPE_TRAITS 1
129		#endif
130		#ifndef HALF_ENABLE_CPP11_CSTDINT
131			#define HALF_ENABLE_CPP11_CSTDINT 1
132		#endif
133		#ifndef HALF_ENABLE_CPP11_CMATH
134			#define HALF_ENABLE_CPP11_CMATH 1
135		#endif
136		#ifndef HALF_ENABLE_CPP11_HASH
137			#define HALF_ENABLE_CPP11_HASH 1
138		#endif
139		#ifndef HALF_ENABLE_CPP11_CFENV
140			#define HALF_ENABLE_CPP11_CFENV 1
141		#endif
142	#endif
143#elif defined(__GLIBCXX__)									// libstdc++
144	#if defined(__GXX_EXPERIMENTAL_CXX0X__) || __cplusplus >= 201103
145		#ifdef __clang__
146			#if __GLIBCXX__ >= 20080606 && !defined(HALF_ENABLE_CPP11_TYPE_TRAITS)
147				#define HALF_ENABLE_CPP11_TYPE_TRAITS 1
148			#endif
149			#if __GLIBCXX__ >= 20080606 && !defined(HALF_ENABLE_CPP11_CSTDINT)
150				#define HALF_ENABLE_CPP11_CSTDINT 1
151			#endif
152			#if __GLIBCXX__ >= 20080606 && !defined(HALF_ENABLE_CPP11_CMATH)
153				#define HALF_ENABLE_CPP11_CMATH 1
154			#endif
155			#if __GLIBCXX__ >= 20080606 && !defined(HALF_ENABLE_CPP11_HASH)
156				#define HALF_ENABLE_CPP11_HASH 1
157			#endif
158			#if __GLIBCXX__ >= 20080606 && !defined(HALF_ENABLE_CPP11_CFENV)
159				#define HALF_ENABLE_CPP11_CFENV 1
160			#endif
161		#else
162			#if HALF_GCC_VERSION >= 403 && !defined(HALF_ENABLE_CPP11_TYPE_TRAITS)
163				#define HALF_ENABLE_CPP11_TYPE_TRAITS 1
164			#endif
165			#if HALF_GCC_VERSION >= 403 && !defined(HALF_ENABLE_CPP11_CSTDINT)
166				#define HALF_ENABLE_CPP11_CSTDINT 1
167			#endif
168			#if HALF_GCC_VERSION >= 403 && !defined(HALF_ENABLE_CPP11_CMATH)
169				#define HALF_ENABLE_CPP11_CMATH 1
170			#endif
171			#if HALF_GCC_VERSION >= 403 && !defined(HALF_ENABLE_CPP11_HASH)
172				#define HALF_ENABLE_CPP11_HASH 1
173			#endif
174			#if HALF_GCC_VERSION >= 403 && !defined(HALF_ENABLE_CPP11_CFENV)
175				#define HALF_ENABLE_CPP11_CFENV 1
176			#endif
177		#endif
178	#endif
179#elif defined(_CPPLIB_VER)									// Dinkumware/Visual C++
180	#if _CPPLIB_VER >= 520 && !defined(HALF_ENABLE_CPP11_TYPE_TRAITS)
181		#define HALF_ENABLE_CPP11_TYPE_TRAITS 1
182	#endif
183	#if _CPPLIB_VER >= 520 && !defined(HALF_ENABLE_CPP11_CSTDINT)
184			#define HALF_ENABLE_CPP11_CSTDINT 1
185	#endif
186	#if _CPPLIB_VER >= 520 && !defined(HALF_ENABLE_CPP11_HASH)
187		#define HALF_ENABLE_CPP11_HASH 1
188	#endif
189	#if _CPPLIB_VER >= 610 && !defined(HALF_ENABLE_CPP11_CMATH)
190		#define HALF_ENABLE_CPP11_CMATH 1
191	#endif
192	#if _CPPLIB_VER >= 610 && !defined(HALF_ENABLE_CPP11_CFENV)
193		#define HALF_ENABLE_CPP11_CFENV 1
194	#endif
195#endif
196#undef HALF_GCC_VERSION
197#undef HALF_ICC_VERSION
198
199// any error throwing C++ exceptions?
200#if defined(HALF_ERRHANDLING_THROW_INVALID) || defined(HALF_ERRHANDLING_THROW_DIVBYZERO) || defined(HALF_ERRHANDLING_THROW_OVERFLOW) || defined(HALF_ERRHANDLING_THROW_UNDERFLOW) || defined(HALF_ERRHANDLING_THROW_INEXACT)
201#define HALF_ERRHANDLING_THROWS 1
202#endif
203
204// any error handling enabled?
205#define HALF_ERRHANDLING	(HALF_ERRHANDLING_FLAGS||HALF_ERRHANDLING_ERRNO||HALF_ERRHANDLING_FENV||HALF_ERRHANDLING_THROWS)
206
207#if HALF_ERRHANDLING
208	#define HALF_UNUSED_NOERR(name) name
209#else
210	#define HALF_UNUSED_NOERR(name)
211#endif
212
213// support constexpr
214#if HALF_ENABLE_CPP11_CONSTEXPR
215	#define HALF_CONSTEXPR				constexpr
216	#define HALF_CONSTEXPR_CONST		constexpr
217	#if HALF_ERRHANDLING
218		#define HALF_CONSTEXPR_NOERR
219	#else
220		#define HALF_CONSTEXPR_NOERR	constexpr
221	#endif
222#else
223	#define HALF_CONSTEXPR
224	#define HALF_CONSTEXPR_CONST		const
225	#define HALF_CONSTEXPR_NOERR
226#endif
227
228// support noexcept
229#if HALF_ENABLE_CPP11_NOEXCEPT
230	#define HALF_NOEXCEPT	noexcept
231	#define HALF_NOTHROW	noexcept
232#else
233	#define HALF_NOEXCEPT
234	#define HALF_NOTHROW	throw()
235#endif
236
237// support thread storage
238#if HALF_ENABLE_CPP11_THREAD_LOCAL
239	#define HALF_THREAD_LOCAL	thread_local
240#else
241	#define HALF_THREAD_LOCAL	static
242#endif
243
244#include <utility>
245#include <algorithm>
246#include <istream>
247#include <ostream>
248#include <limits>
249#include <stdexcept>
250#include <climits>
251#include <cmath>
252#include <cstring>
253#include <cstdlib>
254#if HALF_ENABLE_CPP11_TYPE_TRAITS
255	#include <type_traits>
256#endif
257#if HALF_ENABLE_CPP11_CSTDINT
258	#include <cstdint>
259#endif
260#if HALF_ERRHANDLING_ERRNO
261	#include <cerrno>
262#endif
263#if HALF_ENABLE_CPP11_CFENV
264	#include <cfenv>
265#endif
266#if HALF_ENABLE_CPP11_HASH
267	#include <functional>
268#endif
269
270
271#ifndef HALF_ENABLE_F16C_INTRINSICS
272	/// Enable F16C intruction set intrinsics.
273	/// Defining this to 1 enables the use of [F16C compiler intrinsics](https://en.wikipedia.org/wiki/F16C) for converting between 
274	/// half-precision and single-precision values which may result in improved performance. This will not perform additional checks 
275	/// for support of the F16C instruction set, so an appropriate target platform is required when enabling this feature.
276	///
277	/// Unless predefined it will be enabled automatically when the `__F16C__` symbol is defined, which some compilers do on supporting platforms.
278	#define HALF_ENABLE_F16C_INTRINSICS __F16C__
279#endif
280#if HALF_ENABLE_F16C_INTRINSICS
281	#include <immintrin.h>
282#endif
283
284#ifdef HALF_DOXYGEN_ONLY
285/// Type for internal floating-point computations.
286/// This can be predefined to a built-in floating-point type (`float`, `double` or `long double`) to override the internal 
287/// half-precision implementation to use this type for computing arithmetic operations and mathematical function (if available). 
288/// This can result in improved performance for arithmetic operators and mathematical functions but might cause results to 
289/// deviate from the specified half-precision rounding mode and inhibits proper detection of half-precision exceptions.
290#define HALF_ARITHMETIC_TYPE (undefined)
291
292/// Enable internal exception flags.
293/// Defining this to 1 causes operations on half-precision values to raise internal floating-point exception flags according to 
294/// the IEEE 754 standard. These can then be cleared and checked with clearexcept(), testexcept().
295#define HALF_ERRHANDLING_FLAGS	0
296
297/// Enable exception propagation to `errno`.
298/// Defining this to 1 causes operations on half-precision values to propagate floating-point exceptions to 
299/// [errno](https://en.cppreference.com/w/cpp/error/errno) from `<cerrno>`. Specifically this will propagate domain errors as 
300/// [EDOM](https://en.cppreference.com/w/cpp/error/errno_macros) and pole, overflow and underflow errors as 
301/// [ERANGE](https://en.cppreference.com/w/cpp/error/errno_macros). Inexact errors won't be propagated.
302#define HALF_ERRHANDLING_ERRNO	0
303
304/// Enable exception propagation to built-in floating-point platform.
305/// Defining this to 1 causes operations on half-precision values to propagate floating-point exceptions to the built-in 
306/// single- and double-precision implementation's exception flags using the 
307/// [C++11 floating-point environment control](https://en.cppreference.com/w/cpp/numeric/fenv) from `<cfenv>`. However, this 
308/// does not work in reverse and single- or double-precision exceptions will not raise the corresponding half-precision 
309/// exception flags, nor will explicitly clearing flags clear the corresponding built-in flags.
310#define HALF_ERRHANDLING_FENV	0
311
312/// Throw C++ exception on domain errors.
313/// Defining this to a string literal causes operations on half-precision values to throw a 
314/// [std::domain_error](https://en.cppreference.com/w/cpp/error/domain_error) with the specified message on domain errors.
315#define HALF_ERRHANDLING_THROW_INVALID		(undefined)
316
317/// Throw C++ exception on pole errors.
318/// Defining this to a string literal causes operations on half-precision values to throw a 
319/// [std::domain_error](https://en.cppreference.com/w/cpp/error/domain_error) with the specified message on pole errors.
320#define HALF_ERRHANDLING_THROW_DIVBYZERO	(undefined)
321
322/// Throw C++ exception on overflow errors.
323/// Defining this to a string literal causes operations on half-precision values to throw a 
324/// [std::overflow_error](https://en.cppreference.com/w/cpp/error/overflow_error) with the specified message on overflows.
325#define HALF_ERRHANDLING_THROW_OVERFLOW		(undefined)
326
327/// Throw C++ exception on underflow errors.
328/// Defining this to a string literal causes operations on half-precision values to throw a 
329/// [std::underflow_error](https://en.cppreference.com/w/cpp/error/underflow_error) with the specified message on underflows.
330#define HALF_ERRHANDLING_THROW_UNDERFLOW	(undefined)
331
332/// Throw C++ exception on rounding errors.
333/// Defining this to 1 causes operations on half-precision values to throw a 
334/// [std::range_error](https://en.cppreference.com/w/cpp/error/range_error) with the specified message on general rounding errors.
335#define HALF_ERRHANDLING_THROW_INEXACT		(undefined)
336#endif
337
338#ifndef HALF_ERRHANDLING_OVERFLOW_TO_INEXACT
339/// Raise INEXACT exception on overflow.
340/// Defining this to 1 (default) causes overflow errors to automatically raise inexact exceptions in addition.
341/// These will be raised after any possible handling of the underflow exception.
342#define HALF_ERRHANDLING_OVERFLOW_TO_INEXACT	1
343#endif
344
345#ifndef HALF_ERRHANDLING_UNDERFLOW_TO_INEXACT
346/// Raise INEXACT exception on underflow.
347/// Defining this to 1 (default) causes underflow errors to automatically raise inexact exceptions in addition.
348/// These will be raised after any possible handling of the underflow exception.
349///
350/// **Note:** This will actually cause underflow (and the accompanying inexact) exceptions to be raised *only* when the result 
351/// is inexact, while if disabled bare underflow errors will be raised for *any* (possibly exact) subnormal result.
352#define HALF_ERRHANDLING_UNDERFLOW_TO_INEXACT	1
353#endif
354
355/// Default rounding mode.
356/// This specifies the rounding mode used for all conversions between [half](\ref half_float::half)s and more precise types 
357/// (unless using half_cast() and specifying the rounding mode directly) as well as in arithmetic operations and mathematical 
358/// functions. It can be redefined (before including half.hpp) to one of the standard rounding modes using their respective 
359/// constants or the equivalent values of 
360/// [std::float_round_style](https://en.cppreference.com/w/cpp/types/numeric_limits/float_round_style):
361///
362/// `std::float_round_style`         | value | rounding
363/// ---------------------------------|-------|-------------------------
364/// `std::round_indeterminate`       | -1    | fastest
365/// `std::round_toward_zero`         | 0     | toward zero
366/// `std::round_to_nearest`          | 1     | to nearest (default)
367/// `std::round_toward_infinity`     | 2     | toward positive infinity
368/// `std::round_toward_neg_infinity` | 3     | toward negative infinity
369///
370/// By default this is set to `1` (`std::round_to_nearest`), which rounds results to the nearest representable value. It can even 
371/// be set to [std::numeric_limits<float>::round_style](https://en.cppreference.com/w/cpp/types/numeric_limits/round_style) to synchronize 
372/// the rounding mode with that of the built-in single-precision implementation (which is likely `std::round_to_nearest`, though).
373#ifndef HALF_ROUND_STYLE
374	#define HALF_ROUND_STYLE	1		// = std::round_to_nearest
375#endif
376
377/// Value signaling overflow.
378/// In correspondence with `HUGE_VAL[F|L]` from `<cmath>` this symbol expands to a positive value signaling the overflow of an 
379/// operation, in particular it just evaluates to positive infinity.
380///
381/// **See also:** Documentation for [HUGE_VAL](https://en.cppreference.com/w/cpp/numeric/math/HUGE_VAL)
382#define HUGE_VALH	std::numeric_limits<half_float::half>::infinity()
383
384/// Fast half-precision fma function.
385/// This symbol is defined if the fma() function generally executes as fast as, or faster than, a separate 
386/// half-precision multiplication followed by an addition, which is always the case.
387///
388/// **See also:** Documentation for [FP_FAST_FMA](https://en.cppreference.com/w/cpp/numeric/math/fma)
389#define FP_FAST_FMAH	1
390
391///	Half rounding mode.
392/// In correspondence with `FLT_ROUNDS` from `<cfloat>` this symbol expands to the rounding mode used for 
393/// half-precision operations. It is an alias for [HALF_ROUND_STYLE](\ref HALF_ROUND_STYLE).
394///
395/// **See also:** Documentation for [FLT_ROUNDS](https://en.cppreference.com/w/cpp/types/climits/FLT_ROUNDS)
396#define HLF_ROUNDS	HALF_ROUND_STYLE
397
398#ifndef FP_ILOGB0
399	#define FP_ILOGB0		INT_MIN
400#endif
401#ifndef FP_ILOGBNAN
402	#define FP_ILOGBNAN		INT_MAX
403#endif
404#ifndef FP_SUBNORMAL
405	#define FP_SUBNORMAL	0
406#endif
407#ifndef FP_ZERO
408	#define FP_ZERO			1
409#endif
410#ifndef FP_NAN
411	#define FP_NAN			2
412#endif
413#ifndef FP_INFINITE
414	#define FP_INFINITE		3
415#endif
416#ifndef FP_NORMAL
417	#define FP_NORMAL		4
418#endif
419
420#if !HALF_ENABLE_CPP11_CFENV && !defined(FE_ALL_EXCEPT)
421	#define FE_INVALID		0x10
422	#define FE_DIVBYZERO	0x08
423	#define FE_OVERFLOW		0x04
424	#define FE_UNDERFLOW	0x02
425	#define FE_INEXACT		0x01
426	#define FE_ALL_EXCEPT	(FE_INVALID|FE_DIVBYZERO|FE_OVERFLOW|FE_UNDERFLOW|FE_INEXACT)
427#endif
428
429
430/// Main namespace for half-precision functionality.
431/// This namespace contains all the functionality provided by the library.
432namespace half_float
433{
434	class half;
435
436#if HALF_ENABLE_CPP11_USER_LITERALS
437	/// Library-defined half-precision literals.
438	/// Import this namespace to enable half-precision floating-point literals:
439	/// ~~~~{.cpp}
440	/// using namespace half_float::literal;
441	/// half_float::half = 4.2_h;
442	/// ~~~~
443	namespace literal
444	{
445		half operator "" _h(long double);
446	}
447#endif
448
449	/// \internal
450	/// \brief Implementation details.
451	namespace detail
452	{
453	#if HALF_ENABLE_CPP11_TYPE_TRAITS
454		/// Conditional type.
455		template<bool B,typename T,typename F> struct conditional : std::conditional<B,T,F> {};
456
457		/// Helper for tag dispatching.
458		template<bool B> struct bool_type : std::integral_constant<bool,B> {};
459		using std::true_type;
460		using std::false_type;
461
462		/// Type traits for floating-point types.
463		template<typename T> struct is_float : std::is_floating_point<T> {};
464	#else
465		/// Conditional type.
466		template<bool,typename T,typename> struct conditional { typedef T type; };
467		template<typename T,typename F> struct conditional<false,T,F> { typedef F type; };
468
469		/// Helper for tag dispatching.
470		template<bool> struct bool_type {};
471		typedef bool_type<true> true_type;
472		typedef bool_type<false> false_type;
473
474		/// Type traits for floating-point types.
475		template<typename> struct is_float : false_type {};
476		template<typename T> struct is_float<const T> : is_float<T> {};
477		template<typename T> struct is_float<volatile T> : is_float<T> {};
478		template<typename T> struct is_float<const volatile T> : is_float<T> {};
479		template<> struct is_float<float> : true_type {};
480		template<> struct is_float<double> : true_type {};
481		template<> struct is_float<long double> : true_type {};
482	#endif
483
484		/// Type traits for floating-point bits.
485		template<typename T> struct bits { typedef unsigned char type; };
486		template<typename T> struct bits<const T> : bits<T> {};
487		template<typename T> struct bits<volatile T> : bits<T> {};
488		template<typename T> struct bits<const volatile T> : bits<T> {};
489
490	#if HALF_ENABLE_CPP11_CSTDINT
491		/// Unsigned integer of (at least) 16 bits width.
492		typedef std::uint_least16_t uint16;
493
494		/// Fastest unsigned integer of (at least) 32 bits width.
495		typedef std::uint_fast32_t uint32;
496
497		/// Fastest signed integer of (at least) 32 bits width.
498		typedef std::int_fast32_t int32;
499
500		/// Unsigned integer of (at least) 32 bits width.
501		template<> struct bits<float> { typedef std::uint_least32_t type; };
502
503		/// Unsigned integer of (at least) 64 bits width.
504		template<> struct bits<double> { typedef std::uint_least64_t type; };
505	#else
506		/// Unsigned integer of (at least) 16 bits width.
507		typedef unsigned short uint16;
508
509		/// Fastest unsigned integer of (at least) 32 bits width.
510		typedef unsigned long uint32;
511
512		/// Fastest unsigned integer of (at least) 32 bits width.
513		typedef long int32;
514
515		/// Unsigned integer of (at least) 32 bits width.
516		template<> struct bits<float> : conditional<std::numeric_limits<unsigned int>::digits>=32,unsigned int,unsigned long> {};
517
518		#if HALF_ENABLE_CPP11_LONG_LONG
519			/// Unsigned integer of (at least) 64 bits width.
520			template<> struct bits<double> : conditional<std::numeric_limits<unsigned long>::digits>=64,unsigned long,unsigned long long> {};
521		#else
522			/// Unsigned integer of (at least) 64 bits width.
523			template<> struct bits<double> { typedef unsigned long type; };
524		#endif
525	#endif
526
527	#ifdef HALF_ARITHMETIC_TYPE
528		/// Type to use for arithmetic computations and mathematic functions internally.
529		typedef HALF_ARITHMETIC_TYPE internal_t;
530	#endif
531
532		/// Tag type for binary construction.
533		struct binary_t {};
534
535		/// Tag for binary construction.
536		HALF_CONSTEXPR_CONST binary_t binary = binary_t();
537
538		/// \name Implementation defined classification and arithmetic
539		/// \{
540
541		/// Check for infinity.
542		/// \tparam T argument type (builtin floating-point type)
543		/// \param arg value to query
544		/// \retval true if infinity
545		/// \retval false else
546		template<typename T> bool builtin_isinf(T arg)
547		{
548		#if HALF_ENABLE_CPP11_CMATH
549			return std::isinf(arg);
550		#elif defined(_MSC_VER)
551			return !::_finite(static_cast<double>(arg)) && !::_isnan(static_cast<double>(arg));
552		#else
553			return arg == std::numeric_limits<T>::infinity() || arg == -std::numeric_limits<T>::infinity();
554		#endif
555		}
556
557		/// Check for NaN.
558		/// \tparam T argument type (builtin floating-point type)
559		/// \param arg value to query
560		/// \retval true if not a number
561		/// \retval false else
562		template<typename T> bool builtin_isnan(T arg)
563		{
564		#if HALF_ENABLE_CPP11_CMATH
565			return std::isnan(arg);
566		#elif defined(_MSC_VER)
567			return ::_isnan(static_cast<double>(arg)) != 0;
568		#else
569			return arg != arg;
570		#endif
571		}
572
573		/// Check sign.
574		/// \tparam T argument type (builtin floating-point type)
575		/// \param arg value to query
576		/// \retval true if signbit set
577		/// \retval false else
578		template<typename T> bool builtin_signbit(T arg)
579		{
580		#if HALF_ENABLE_CPP11_CMATH
581			return std::signbit(arg);
582		#else
583			return arg < T() || (arg == T() && T(1)/arg < T());
584		#endif
585		}
586
587		/// Platform-independent sign mask.
588		/// \param arg integer value in two's complement
589		/// \retval -1 if \a arg negative
590		/// \retval 0 if \a arg positive
591		inline uint32 sign_mask(uint32 arg)
592		{
593			static const int N = std::numeric_limits<uint32>::digits - 1;
594		#if HALF_TWOS_COMPLEMENT_INT
595			return static_cast<int32>(arg) >> N;
596		#else
597			return -((arg>>N)&1);
598		#endif
599		}
600
601		/// Platform-independent arithmetic right shift.
602		/// \param arg integer value in two's complement
603		/// \param i shift amount (at most 31)
604		/// \return \a arg right shifted for \a i bits with possible sign extension
605		inline uint32 arithmetic_shift(uint32 arg, int i)
606		{
607		#if HALF_TWOS_COMPLEMENT_INT
608			return static_cast<int32>(arg) >> i;
609		#else
610			return static_cast<int32>(arg)/(static_cast<int32>(1)<<i) - ((arg>>(std::numeric_limits<uint32>::digits-1))&1);
611		#endif
612		}
613
614		/// \}
615		/// \name Error handling
616		/// \{
617
618		/// Internal exception flags.
619		/// \return reference to global exception flags
620		inline int& errflags() { HALF_THREAD_LOCAL int flags = 0; return flags; }
621
622		/// Raise floating-point exception.
623		/// \param flags exceptions to raise
624		/// \param cond condition to raise exceptions for
625		inline void raise(int HALF_UNUSED_NOERR(flags), bool HALF_UNUSED_NOERR(cond) = true)
626		{
627		#if HALF_ERRHANDLING
628			if(!cond)
629				return;
630		#if HALF_ERRHANDLING_FLAGS
631			errflags() |= flags;
632		#endif
633		#if HALF_ERRHANDLING_ERRNO
634			if(flags & FE_INVALID)
635				errno = EDOM;
636			else if(flags & (FE_DIVBYZERO|FE_OVERFLOW|FE_UNDERFLOW))
637				errno = ERANGE;
638		#endif
639		#if HALF_ERRHANDLING_FENV && HALF_ENABLE_CPP11_CFENV
640			std::feraiseexcept(flags);
641		#endif
642		#ifdef HALF_ERRHANDLING_THROW_INVALID
643			if(flags & FE_INVALID)
644				throw std::domain_error(HALF_ERRHANDLING_THROW_INVALID);
645		#endif
646		#ifdef HALF_ERRHANDLING_THROW_DIVBYZERO
647			if(flags & FE_DIVBYZERO)
648				throw std::domain_error(HALF_ERRHANDLING_THROW_DIVBYZERO);
649		#endif
650		#ifdef HALF_ERRHANDLING_THROW_OVERFLOW
651			if(flags & FE_OVERFLOW)
652				throw std::overflow_error(HALF_ERRHANDLING_THROW_OVERFLOW);
653		#endif
654		#ifdef HALF_ERRHANDLING_THROW_UNDERFLOW
655			if(flags & FE_UNDERFLOW)
656				throw std::underflow_error(HALF_ERRHANDLING_THROW_UNDERFLOW);
657		#endif
658		#ifdef HALF_ERRHANDLING_THROW_INEXACT
659			if(flags & FE_INEXACT)
660				throw std::range_error(HALF_ERRHANDLING_THROW_INEXACT);
661		#endif
662		#if HALF_ERRHANDLING_UNDERFLOW_TO_INEXACT
663			if((flags & FE_UNDERFLOW) && !(flags & FE_INEXACT))
664				raise(FE_INEXACT);
665		#endif
666		#if HALF_ERRHANDLING_OVERFLOW_TO_INEXACT
667			if((flags & FE_OVERFLOW) && !(flags & FE_INEXACT))
668				raise(FE_INEXACT);
669		#endif
670		#endif
671		}
672
673		/// Check and signal for any NaN.
674		/// \param x first half-precision value to check
675		/// \param y second half-precision value to check
676		/// \retval true if either \a x or \a y is NaN
677		/// \retval false else
678		/// \exception FE_INVALID if \a x or \a y is NaN
679		inline HALF_CONSTEXPR_NOERR bool compsignal(unsigned int x, unsigned int y)
680		{
681		#if HALF_ERRHANDLING
682			raise(FE_INVALID, (x&0x7FFF)>0x7C00 || (y&0x7FFF)>0x7C00);
683		#endif
684			return (x&0x7FFF) > 0x7C00 || (y&0x7FFF) > 0x7C00;
685		}
686
687		/// Signal and silence signaling NaN.
688		/// \param nan half-precision NaN value
689		/// \return quiet NaN
690		/// \exception FE_INVALID if \a nan is signaling NaN
691		inline HALF_CONSTEXPR_NOERR unsigned int signal(unsigned int nan)
692		{
693		#if HALF_ERRHANDLING
694			raise(FE_INVALID, !(nan&0x200));
695		#endif
696			return nan | 0x200;
697		}
698
699		/// Signal and silence signaling NaNs.
700		/// \param x first half-precision value to check
701		/// \param y second half-precision value to check
702		/// \return quiet NaN
703		/// \exception FE_INVALID if \a x or \a y is signaling NaN
704		inline HALF_CONSTEXPR_NOERR unsigned int signal(unsigned int x, unsigned int y)
705		{
706		#if HALF_ERRHANDLING
707			raise(FE_INVALID, ((x&0x7FFF)>0x7C00 && !(x&0x200)) || ((y&0x7FFF)>0x7C00 && !(y&0x200)));
708		#endif
709			return ((x&0x7FFF)>0x7C00) ? (x|0x200) : (y|0x200);
710		}
711
712		/// Signal and silence signaling NaNs.
713		/// \param x first half-precision value to check
714		/// \param y second half-precision value to check
715		/// \param z third half-precision value to check
716		/// \return quiet NaN
717		/// \exception FE_INVALID if \a x, \a y or \a z is signaling NaN
718		inline HALF_CONSTEXPR_NOERR unsigned int signal(unsigned int x, unsigned int y, unsigned int z)
719		{
720		#if HALF_ERRHANDLING
721			raise(FE_INVALID, ((x&0x7FFF)>0x7C00 && !(x&0x200)) || ((y&0x7FFF)>0x7C00 && !(y&0x200)) || ((z&0x7FFF)>0x7C00 && !(z&0x200)));
722		#endif
723			return ((x&0x7FFF)>0x7C00) ? (x|0x200) : ((y&0x7FFF)>0x7C00) ? (y|0x200) : (z|0x200);
724		}
725
726		/// Select value or signaling NaN.
727		/// \param x preferred half-precision value
728		/// \param y ignored half-precision value except for signaling NaN
729		/// \return \a y if signaling NaN, \a x otherwise
730		/// \exception FE_INVALID if \a y is signaling NaN
731		inline HALF_CONSTEXPR_NOERR unsigned int select(unsigned int x, unsigned int HALF_UNUSED_NOERR(y))
732		{
733		#if HALF_ERRHANDLING
734			return (((y&0x7FFF)>0x7C00) && !(y&0x200)) ? signal(y) : x;
735		#else
736			return x;
737		#endif
738		}
739
740		/// Raise domain error and return NaN.
741		/// return quiet NaN
742		/// \exception FE_INVALID
743		inline HALF_CONSTEXPR_NOERR unsigned int invalid()
744		{
745		#if HALF_ERRHANDLING
746			raise(FE_INVALID);
747		#endif
748			return 0x7FFF;
749		}
750
751		/// Raise pole error and return infinity.
752		/// \param sign half-precision value with sign bit only
753		/// \return half-precision infinity with sign of \a sign
754		/// \exception FE_DIVBYZERO
755		inline HALF_CONSTEXPR_NOERR unsigned int pole(unsigned int sign = 0)
756		{
757		#if HALF_ERRHANDLING
758			raise(FE_DIVBYZERO);
759		#endif
760			return sign | 0x7C00;
761		}
762
763		/// Check value for underflow.
764		/// \param arg non-zero half-precision value to check
765		/// \return \a arg
766		/// \exception FE_UNDERFLOW if arg is subnormal
767		inline HALF_CONSTEXPR_NOERR unsigned int check_underflow(unsigned int arg)
768		{
769		#if HALF_ERRHANDLING && !HALF_ERRHANDLING_UNDERFLOW_TO_INEXACT
770			raise(FE_UNDERFLOW, !(arg&0x7C00));
771		#endif
772			return arg;
773		}
774
775		/// \}
776		/// \name Conversion and rounding
777		/// \{
778
779		/// Half-precision overflow.
780		/// \tparam R rounding mode to use
781		/// \param sign half-precision value with sign bit only
782		/// \return rounded overflowing half-precision value
783		/// \exception FE_OVERFLOW
784		template<std::float_round_style R> HALF_CONSTEXPR_NOERR unsigned int overflow(unsigned int sign = 0)
785		{
786		#if HALF_ERRHANDLING
787			raise(FE_OVERFLOW);
788		#endif
789			return	(R==std::round_toward_infinity) ? (sign+0x7C00-(sign>>15)) :
790					(R==std::round_toward_neg_infinity) ? (sign+0x7BFF+(sign>>15)) :
791					(R==std::round_toward_zero) ? (sign|0x7BFF) :
792					(sign|0x7C00);
793		}
794
795		/// Half-precision underflow.
796		/// \tparam R rounding mode to use
797		/// \param sign half-precision value with sign bit only
798		/// \return rounded underflowing half-precision value
799		/// \exception FE_UNDERFLOW
800		template<std::float_round_style R> HALF_CONSTEXPR_NOERR unsigned int underflow(unsigned int sign = 0)
801		{
802		#if HALF_ERRHANDLING
803			raise(FE_UNDERFLOW);
804		#endif
805			return	(R==std::round_toward_infinity) ? (sign+1-(sign>>15)) :
806					(R==std::round_toward_neg_infinity) ? (sign+(sign>>15)) :
807					sign;
808		}
809
810		/// Round half-precision number.
811		/// \tparam R rounding mode to use
812		/// \tparam I `true` to always raise INEXACT exception, `false` to raise only for rounded results
813		/// \param value finite half-precision number to round
814		/// \param g guard bit (most significant discarded bit)
815		/// \param s sticky bit (or of all but the most significant discarded bits)
816		/// \return rounded half-precision value
817		/// \exception FE_OVERFLOW on overflows
818		/// \exception FE_UNDERFLOW on underflows
819		/// \exception FE_INEXACT if value had to be rounded or \a I is `true`
820		template<std::float_round_style R,bool I> HALF_CONSTEXPR_NOERR unsigned int rounded(unsigned int value, int g, int s)
821		{
822		#if HALF_ERRHANDLING
823			value +=	(R==std::round_to_nearest) ? (g&(s|value)) :
824						(R==std::round_toward_infinity) ? (~(value>>15)&(g|s)) :
825						(R==std::round_toward_neg_infinity) ? ((value>>15)&(g|s)) : 0;
826			if((value&0x7C00) == 0x7C00)
827				raise(FE_OVERFLOW);
828			else if(value & 0x7C00)
829				raise(FE_INEXACT, I || (g|s)!=0);
830			else
831				raise(FE_UNDERFLOW, !(HALF_ERRHANDLING_UNDERFLOW_TO_INEXACT) || I || (g|s)!=0);
832			return value;
833		#else
834			return	(R==std::round_to_nearest) ? (value+(g&(s|value))) :
835					(R==std::round_toward_infinity) ? (value+(~(value>>15)&(g|s))) :
836					(R==std::round_toward_neg_infinity) ? (value+((value>>15)&(g|s))) :
837					value;
838		#endif
839		}
840
841		/// Round half-precision number to nearest integer value.
842		/// \tparam R rounding mode to use
843		/// \tparam E `true` for round to even, `false` for round away from zero
844		/// \tparam I `true` to raise INEXACT exception (if inexact), `false` to never raise it
845		/// \param value half-precision value to round
846		/// \return half-precision bits for nearest integral value
847		/// \exception FE_INVALID for signaling NaN
848		/// \exception FE_INEXACT if value had to be rounded and \a I is `true`
849		template<std::float_round_style R,bool E,bool I> unsigned int integral(unsigned int value)
850		{
851			unsigned int abs = value & 0x7FFF;
852			if(abs < 0x3C00)
853			{
854				raise(FE_INEXACT, I);
855				return ((R==std::round_to_nearest) ? (0x3C00&-static_cast<unsigned>(abs>=(0x3800+E))) :
856						(R==std::round_toward_infinity) ? (0x3C00&-(~(value>>15)&(abs!=0))) :
857						(R==std::round_toward_neg_infinity) ? (0x3C00&-static_cast<unsigned>(value>0x8000)) :
858						0) | (value&0x8000);
859			}
860			if(abs >= 0x6400)
861				return (abs>0x7C00) ? signal(value) : value;
862			unsigned int exp = 25 - (abs>>10), mask = (1<<exp) - 1;
863			raise(FE_INEXACT, I && (value&mask));
864			return ((	(R==std::round_to_nearest) ? ((1<<(exp-1))-(~(value>>exp)&E)) :
865						(R==std::round_toward_infinity) ? (mask&((value>>15)-1)) :
866						(R==std::round_toward_neg_infinity) ? (mask&-(value>>15)) :
867						0) + value) & ~mask;
868		}
869
870		/// Convert fixed point to half-precision floating-point.
871		/// \tparam R rounding mode to use
872		/// \tparam F number of fractional bits in [11,31]
873		/// \tparam S `true` for signed, `false` for unsigned
874		/// \tparam N `true` for additional normalization step, `false` if already normalized to 1.F
875		/// \tparam I `true` to always raise INEXACT exception, `false` to raise only for rounded results
876		/// \param m mantissa in Q1.F fixed point format
877		/// \param exp biased exponent - 1
878		/// \param sign half-precision value with sign bit only
879		/// \param s sticky bit (or of all but the most significant already discarded bits)
880		/// \return value converted to half-precision
881		/// \exception FE_OVERFLOW on overflows
882		/// \exception FE_UNDERFLOW on underflows
883		/// \exception FE_INEXACT if value had to be rounded or \a I is `true`
884		template<std::float_round_style R,unsigned int F,bool S,bool N,bool I> unsigned int fixed2half(uint32 m, int exp = 14, unsigned int sign = 0, int s = 0)
885		{
886			if(S)
887			{
888				uint32 msign = sign_mask(m);
889				m = (m^msign) - msign;
890				sign = msign & 0x8000;
891			}
892			if(N)
893				for(; m<(static_cast<uint32>(1)<<F) && exp; m<<=1,--exp) ;
894			else if(exp < 0)
895				return rounded<R,I>(sign+(m>>(F-10-exp)), (m>>(F-11-exp))&1, s|((m&((static_cast<uint32>(1)<<(F-11-exp))-1))!=0));
896			return rounded<R,I>(sign+(exp<<10)+(m>>(F-10)), (m>>(F-11))&1, s|((m&((static_cast<uint32>(1)<<(F-11))-1))!=0));
897		}
898
899		/// Convert IEEE single-precision to half-precision.
900		/// Credit for this goes to [Jeroen van der Zijp](ftp://ftp.fox-toolkit.org/pub/fasthalffloatconversion.pdf).
901		/// \tparam R rounding mode to use
902		/// \param value single-precision value to convert
903		/// \return rounded half-precision value
904		/// \exception FE_OVERFLOW on overflows
905		/// \exception FE_UNDERFLOW on underflows
906		/// \exception FE_INEXACT if value had to be rounded
907		template<std::float_round_style R> unsigned int float2half_impl(float value, true_type)
908		{
909		#if HALF_ENABLE_F16C_INTRINSICS
910			return _mm_cvtsi128_si32(_mm_cvtps_ph(_mm_set_ss(value),
911				(R==std::round_to_nearest) ? _MM_FROUND_TO_NEAREST_INT :
912				(R==std::round_toward_zero) ? _MM_FROUND_TO_ZERO :
913				(R==std::round_toward_infinity) ? _MM_FROUND_TO_POS_INF :
914				(R==std::round_toward_neg_infinity) ? _MM_FROUND_TO_NEG_INF :
915				_MM_FROUND_CUR_DIRECTION));
916		#else
917			bits<float>::type fbits;
918			std::memcpy(&fbits, &value, sizeof(float));
919		#if 1
920			unsigned int sign = (fbits>>16) & 0x8000;
921			fbits &= 0x7FFFFFFF;
922			if(fbits >= 0x7F800000)
923				return sign | 0x7C00 | ((fbits>0x7F800000) ? (0x200|((fbits>>13)&0x3FF)) : 0);
924			if(fbits >= 0x47800000)
925				return overflow<R>(sign);
926			if(fbits >= 0x38800000)
927				return rounded<R,false>(sign|(((fbits>>23)-112)<<10)|((fbits>>13)&0x3FF), (fbits>>12)&1, (fbits&0xFFF)!=0);
928			if(fbits >= 0x33000000)
929			{
930				int i = 125 - (fbits>>23);
931				fbits = (fbits&0x7FFFFF) | 0x800000;
932				return rounded<R,false>(sign|(fbits>>(i+1)), (fbits>>i)&1, (fbits&((static_cast<uint32>(1)<<i)-1))!=0);
933			}
934			if(fbits != 0)
935				return underflow<R>(sign);
936			return sign;
937		#else
938			static const uint16 base_table[512] = {
939				0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 
940				0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 
941				0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 
942				0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 
943				0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 
944				0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 
945				0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0001, 0x0002, 0x0004, 0x0008, 0x0010, 0x0020, 0x0040, 0x0080, 0x0100, 
946				0x0200, 0x0400, 0x0800, 0x0C00, 0x1000, 0x1400, 0x1800, 0x1C00, 0x2000, 0x2400, 0x2800, 0x2C00, 0x3000, 0x3400, 0x3800, 0x3C00, 
947				0x4000, 0x4400, 0x4800, 0x4C00, 0x5000, 0x5400, 0x5800, 0x5C00, 0x6000, 0x6400, 0x6800, 0x6C00, 0x7000, 0x7400, 0x7800, 0x7BFF, 
948				0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 
949				0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 
950				0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 
951				0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 
952				0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 
953				0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 
954				0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7BFF, 0x7C00, 
955				0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 
956				0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 
957				0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 
958				0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 
959				0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 
960				0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 
961				0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8000, 0x8001, 0x8002, 0x8004, 0x8008, 0x8010, 0x8020, 0x8040, 0x8080, 0x8100, 
962				0x8200, 0x8400, 0x8800, 0x8C00, 0x9000, 0x9400, 0x9800, 0x9C00, 0xA000, 0xA400, 0xA800, 0xAC00, 0xB000, 0xB400, 0xB800, 0xBC00, 
963				0xC000, 0xC400, 0xC800, 0xCC00, 0xD000, 0xD400, 0xD800, 0xDC00, 0xE000, 0xE400, 0xE800, 0xEC00, 0xF000, 0xF400, 0xF800, 0xFBFF, 
964				0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 
965				0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 
966				0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 
967				0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 
968				0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 
969				0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 
970				0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFBFF, 0xFC00 };
971			static const unsigned char shift_table[256] = {
972				24, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 
973				25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 
974				25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 
975				25, 25, 25, 25, 25, 25, 24, 23, 22, 21, 20, 19, 18, 17, 16, 15, 14, 13, 13, 13, 13, 13, 13, 13, 13, 13, 13, 13, 13, 13, 13, 13, 
976				13, 13, 13, 13, 13, 13, 13, 13, 13, 13, 13, 13, 13, 13, 13, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 
977				24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 
978				24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 
979				24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 24, 13 };
980			int sexp = fbits >> 23, exp = sexp & 0xFF, i = shift_table[exp];
981			fbits &= 0x7FFFFF;
982			uint32 m = (fbits|((exp!=0)<<23)) & -static_cast<uint32>(exp!=0xFF);
983			return rounded<R,false>(base_table[sexp]+(fbits>>i), (m>>(i-1))&1, (((static_cast<uint32>(1)<<(i-1))-1)&m)!=0);
984		#endif
985		#endif
986		}
987
988		/// Convert IEEE double-precision to half-precision.
989		/// \tparam R rounding mode to use
990		/// \param value double-precision value to convert
991		/// \return rounded half-precision value
992		/// \exception FE_OVERFLOW on overflows
993		/// \exception FE_UNDERFLOW on underflows
994		/// \exception FE_INEXACT if value had to be rounded
995		template<std::float_round_style R> unsigned int float2half_impl(double value, true_type)
996		{
997		#if HALF_ENABLE_F16C_INTRINSICS
998			if(R == std::round_indeterminate)
999				return _mm_cvtsi128_si32(_mm_cvtps_ph(_mm_cvtpd_ps(_mm_set_sd(value)), _MM_FROUND_CUR_DIRECTION));
1000		#endif
1001			bits<double>::type dbits;
1002			std::memcpy(&dbits, &value, sizeof(double));
1003			uint32 hi = dbits >> 32, lo = dbits & 0xFFFFFFFF;
1004			unsigned int sign = (hi>>16) & 0x8000;
1005			hi &= 0x7FFFFFFF;
1006			if(hi >= 0x7FF00000)
1007				return sign | 0x7C00 | ((dbits&0xFFFFFFFFFFFFF) ? (0x200|((hi>>10)&0x3FF)) : 0);
1008			if(hi >= 0x40F00000)
1009				return overflow<R>(sign);
1010			if(hi >= 0x3F100000)
1011				return rounded<R,false>(sign|(((hi>>20)-1008)<<10)|((hi>>10)&0x3FF), (hi>>9)&1, ((hi&0x1FF)|lo)!=0);
1012			if(hi >= 0x3E600000)
1013			{
1014				int i = 1018 - (hi>>20);
1015				hi = (hi&0xFFFFF) | 0x100000;
1016				return rounded<R,false>(sign|(hi>>(i+1)), (hi>>i)&1, ((hi&((static_cast<uint32>(1)<<i)-1))|lo)!=0);
1017			}
1018			if((hi|lo) != 0)
1019				return underflow<R>(sign);
1020			return sign;
1021		}
1022
1023		/// Convert non-IEEE floating-point to half-precision.
1024		/// \tparam R rounding mode to use
1025		/// \tparam T source type (builtin floating-point type)
1026		/// \param value floating-point value to convert
1027		/// \return rounded half-precision value
1028		/// \exception FE_OVERFLOW on overflows
1029		/// \exception FE_UNDERFLOW on underflows
1030		/// \exception FE_INEXACT if value had to be rounded
1031		template<std::float_round_style R,typename T> unsigned int float2half_impl(T value, ...)
1032		{
1033			unsigned int hbits = static_cast<unsigned>(builtin_signbit(value)) << 15;
1034			if(value == T())
1035				return hbits;
1036			if(builtin_isnan(value))
1037				return hbits | 0x7FFF;
1038			if(builtin_isinf(value))
1039				return hbits | 0x7C00;
1040			int exp;
1041			std::frexp(value, &exp);
1042			if(exp > 16)
1043				return overflow<R>(hbits);
1044			if(exp < -13)
1045				value = std::ldexp(value, 25);
1046			else
1047			{
1048				value = std::ldexp(value, 12-exp);
1049				hbits |= ((exp+13)<<10);
1050			}
1051			T ival, frac = std::modf(value, &ival);
1052			int m = std::abs(static_cast<int>(ival));
1053			return rounded<R,false>(hbits+(m>>1), m&1, frac!=T());
1054		}
1055
1056		/// Convert floating-point to half-precision.
1057		/// \tparam R rounding mode to use
1058		/// \tparam T source type (builtin floating-point type)
1059		/// \param value floating-point value to convert
1060		/// \return rounded half-precision value
1061		/// \exception FE_OVERFLOW on overflows
1062		/// \exception FE_UNDERFLOW on underflows
1063		/// \exception FE_INEXACT if value had to be rounded
1064		template<std::float_round_style R,typename T> unsigned int float2half(T value)
1065		{
1066			return float2half_impl<R>(value, bool_type<std::numeric_limits<T>::is_iec559&&sizeof(typename bits<T>::type)==sizeof(T)>());
1067		}
1068
1069		/// Convert integer to half-precision floating-point.
1070		/// \tparam R rounding mode to use
1071		/// \tparam T type to convert (builtin integer type)
1072		/// \param value integral value to convert
1073		/// \return rounded half-precision value
1074		/// \exception FE_OVERFLOW on overflows
1075		/// \exception FE_INEXACT if value had to be rounded
1076		template<std::float_round_style R,typename T> unsigned int int2half(T value)
1077		{
1078			unsigned int bits = static_cast<unsigned>(value<0) << 15;
1079			if(!value)
1080				return bits;
1081			if(bits)
1082				value = -value;
1083			if(value > 0xFFFF)
1084				return overflow<R>(bits);
1085			unsigned int m = static_cast<unsigned int>(value), exp = 24;
1086			for(; m<0x400; m<<=1,--exp) ;
1087			for(; m>0x7FF; m>>=1,++exp) ;
1088			bits |= (exp<<10) + m;
1089			return (exp>24) ? rounded<R,false>(bits, (value>>(exp-25))&1, (((1<<(exp-25))-1)&value)!=0) : bits;
1090		}
1091
1092		/// Convert half-precision to IEEE single-precision.
1093		/// Credit for this goes to [Jeroen van der Zijp](ftp://ftp.fox-toolkit.org/pub/fasthalffloatconversion.pdf).
1094		/// \param value half-precision value to convert
1095		/// \return single-precision value
1096		inline float half2float_impl(unsigned int value, float, true_type)
1097		{
1098		#if HALF_ENABLE_F16C_INTRINSICS
1099			return _mm_cvtss_f32(_mm_cvtph_ps(_mm_cvtsi32_si128(value)));
1100		#else
1101		#if 0
1102			bits<float>::type fbits = static_cast<bits<float>::type>(value&0x8000) << 16;
1103			int abs = value & 0x7FFF;
1104			if(abs)
1105			{
1106				fbits |= 0x38000000 << static_cast<unsigned>(abs>=0x7C00);
1107				for(; abs<0x400; abs<<=1,fbits-=0x800000) ;
1108				fbits += static_cast<bits<float>::type>(abs) << 13;
1109			}
1110		#else
1111			static const bits<float>::type mantissa_table[2048] = {
1112				0x00000000, 0x33800000, 0x34000000, 0x34400000, 0x34800000, 0x34A00000, 0x34C00000, 0x34E00000, 0x35000000, 0x35100000, 0x35200000, 0x35300000, 0x35400000, 0x35500000, 0x35600000, 0x35700000, 
1113				0x35800000, 0x35880000, 0x35900000, 0x35980000, 0x35A00000, 0x35A80000, 0x35B00000, 0x35B80000, 0x35C00000, 0x35C80000, 0x35D00000, 0x35D80000, 0x35E00000, 0x35E80000, 0x35F00000, 0x35F80000, 
1114				0x36000000, 0x36040000, 0x36080000, 0x360C0000, 0x36100000, 0x36140000, 0x36180000, 0x361C0000, 0x36200000, 0x36240000, 0x36280000, 0x362C0000, 0x36300000, 0x36340000, 0x36380000, 0x363C0000, 
1115				0x36400000, 0x36440000, 0x36480000, 0x364C0000, 0x36500000, 0x36540000, 0x36580000, 0x365C0000, 0x36600000, 0x36640000, 0x36680000, 0x366C0000, 0x36700000, 0x36740000, 0x36780000, 0x367C0000, 
1116				0x36800000, 0x36820000, 0x36840000, 0x36860000, 0x36880000, 0x368A0000, 0x368C0000, 0x368E0000, 0x36900000, 0x36920000, 0x36940000, 0x36960000, 0x36980000, 0x369A0000, 0x369C0000, 0x369E0000, 
1117				0x36A00000, 0x36A20000, 0x36A40000, 0x36A60000, 0x36A80000, 0x36AA0000, 0x36AC0000, 0x36AE0000, 0x36B00000, 0x36B20000, 0x36B40000, 0x36B60000, 0x36B80000, 0x36BA0000, 0x36BC0000, 0x36BE0000, 
1118				0x36C00000, 0x36C20000, 0x36C40000, 0x36C60000, 0x36C80000, 0x36CA0000, 0x36CC0000, 0x36CE0000, 0x36D00000, 0x36D20000, 0x36D40000, 0x36D60000, 0x36D80000, 0x36DA0000, 0x36DC0000, 0x36DE0000, 
1119				0x36E00000, 0x36E20000, 0x36E40000, 0x36E60000, 0x36E80000, 0x36EA0000, 0x36EC0000, 0x36EE0000, 0x36F00000, 0x36F20000, 0x36F40000, 0x36F60000, 0x36F80000, 0x36FA0000, 0x36FC0000, 0x36FE0000, 
1120				0x37000000, 0x37010000, 0x37020000, 0x37030000, 0x37040000, 0x37050000, 0x37060000, 0x37070000, 0x37080000, 0x37090000, 0x370A0000, 0x370B0000, 0x370C0000, 0x370D0000, 0x370E0000, 0x370F0000, 
1121				0x37100000, 0x37110000, 0x37120000, 0x37130000, 0x37140000, 0x37150000, 0x37160000, 0x37170000, 0x37180000, 0x37190000, 0x371A0000, 0x371B0000, 0x371C0000, 0x371D0000, 0x371E0000, 0x371F0000, 
1122				0x37200000, 0x37210000, 0x37220000, 0x37230000, 0x37240000, 0x37250000, 0x37260000, 0x37270000, 0x37280000, 0x37290000, 0x372A0000, 0x372B0000, 0x372C0000, 0x372D0000, 0x372E0000, 0x372F0000, 
1123				0x37300000, 0x37310000, 0x37320000, 0x37330000, 0x37340000, 0x37350000, 0x37360000, 0x37370000, 0x37380000, 0x37390000, 0x373A0000, 0x373B0000, 0x373C0000, 0x373D0000, 0x373E0000, 0x373F0000, 
1124				0x37400000, 0x37410000, 0x37420000, 0x37430000, 0x37440000, 0x37450000, 0x37460000, 0x37470000, 0x37480000, 0x37490000, 0x374A0000, 0x374B0000, 0x374C0000, 0x374D0000, 0x374E0000, 0x374F0000, 
1125				0x37500000, 0x37510000, 0x37520000, 0x37530000, 0x37540000, 0x37550000, 0x37560000, 0x37570000, 0x37580000, 0x37590000, 0x375A0000, 0x375B0000, 0x375C0000, 0x375D0000, 0x375E0000, 0x375F0000, 
1126				0x37600000, 0x37610000, 0x37620000, 0x37630000, 0x37640000, 0x37650000, 0x37660000, 0x37670000, 0x37680000, 0x37690000, 0x376A0000, 0x376B0000, 0x376C0000, 0x376D0000, 0x376E0000, 0x376F0000, 
1127				0x37700000, 0x37710000, 0x37720000, 0x37730000, 0x37740000, 0x37750000, 0x37760000, 0x37770000, 0x37780000, 0x37790000, 0x377A0000, 0x377B0000, 0x377C0000, 0x377D0000, 0x377E0000, 0x377F0000, 
1128				0x37800000, 0x37808000, 0x37810000, 0x37818000, 0x37820000, 0x37828000, 0x37830000, 0x37838000, 0x37840000, 0x37848000, 0x37850000, 0x37858000, 0x37860000, 0x37868000, 0x37870000, 0x37878000, 
1129				0x37880000, 0x37888000, 0x37890000, 0x37898000, 0x378A0000, 0x378A8000, 0x378B0000, 0x378B8000, 0x378C0000, 0x378C8000, 0x378D0000, 0x378D8000, 0x378E0000, 0x378E8000, 0x378F0000, 0x378F8000, 
1130				0x37900000, 0x37908000, 0x37910000, 0x37918000, 0x37920000, 0x37928000, 0x37930000, 0x37938000, 0x37940000, 0x37948000, 0x37950000, 0x37958000, 0x37960000, 0x37968000, 0x37970000, 0x37978000, 
1131				0x37980000, 0x37988000, 0x37990000, 0x37998000, 0x379A0000, 0x379A8000, 0x379B0000, 0x379B8000, 0x379C0000, 0x379C8000, 0x379D0000, 0x379D8000, 0x379E0000, 0x379E8000, 0x379F0000, 0x379F8000, 
1132				0x37A00000, 0x37A08000, 0x37A10000, 0x37A18000, 0x37A20000, 0x37A28000, 0x37A30000, 0x37A38000, 0x37A40000, 0x37A48000, 0x37A50000, 0x37A58000, 0x37A60000, 0x37A68000, 0x37A70000, 0x37A78000, 
1133				0x37A80000, 0x37A88000, 0x37A90000, 0x37A98000, 0x37AA0000, 0x37AA8000, 0x37AB0000, 0x37AB8000, 0x37AC0000, 0x37AC8000, 0x37AD0000, 0x37AD8000, 0x37AE0000, 0x37AE8000, 0x37AF0000, 0x37AF8000, 
1134				0x37B00000, 0x37B08000, 0x37B10000, 0x37B18000, 0x37B20000, 0x37B28000, 0x37B30000, 0x37B38000, 0x37B40000, 0x37B48000, 0x37B50000, 0x37B58000, 0x37B60000, 0x37B68000, 0x37B70000, 0x37B78000, 
1135				0x37B80000, 0x37B88000, 0x37B90000, 0x37B98000, 0x37BA0000, 0x37BA8000, 0x37BB0000, 0x37BB8000, 0x37BC0000, 0x37BC8000, 0x37BD0000, 0x37BD8000, 0x37BE0000, 0x37BE8000, 0x37BF0000, 0x37BF8000, 
1136				0x37C00000, 0x37C08000, 0x37C10000, 0x37C18000, 0x37C20000, 0x37C28000, 0x37C30000, 0x37C38000, 0x37C40000, 0x37C48000, 0x37C50000, 0x37C58000, 0x37C60000, 0x37C68000, 0x37C70000, 0x37C78000, 
1137				0x37C80000, 0x37C88000, 0x37C90000, 0x37C98000, 0x37CA0000, 0x37CA8000, 0x37CB0000, 0x37CB8000, 0x37CC0000, 0x37CC8000, 0x37CD0000, 0x37CD8000, 0x37CE0000, 0x37CE8000, 0x37CF0000, 0x37CF8000, 
1138				0x37D00000, 0x37D08000, 0x37D10000, 0x37D18000, 0x37D20000, 0x37D28000, 0x37D30000, 0x37D38000, 0x37D40000, 0x37D48000, 0x37D50000, 0x37D58000, 0x37D60000, 0x37D68000, 0x37D70000, 0x37D78000, 
1139				0x37D80000, 0x37D88000, 0x37D90000, 0x37D98000, 0x37DA0000, 0x37DA8000, 0x37DB0000, 0x37DB8000, 0x37DC0000, 0x37DC8000, 0x37DD0000, 0x37DD8000, 0x37DE0000, 0x37DE8000, 0x37DF0000, 0x37DF8000, 
1140				0x37E00000, 0x37E08000, 0x37E10000, 0x37E18000, 0x37E20000, 0x37E28000, 0x37E30000, 0x37E38000, 0x37E40000, 0x37E48000, 0x37E50000, 0x37E58000, 0x37E60000, 0x37E68000, 0x37E70000, 0x37E78000, 
1141				0x37E80000, 0x37E88000, 0x37E90000, 0x37E98000, 0x37EA0000, 0x37EA8000, 0x37EB0000, 0x37EB8000, 0x37EC0000, 0x37EC8000, 0x37ED0000, 0x37ED8000, 0x37EE0000, 0x37EE8000, 0x37EF0000, 0x37EF8000, 
1142				0x37F00000, 0x37F08000, 0x37F10000, 0x37F18000, 0x37F20000, 0x37F28000, 0x37F30000, 0x37F38000, 0x37F40000, 0x37F48000, 0x37F50000, 0x37F58000, 0x37F60000, 0x37F68000, 0x37F70000, 0x37F78000, 
1143				0x37F80000, 0x37F88000, 0x37F90000, 0x37F98000, 0x37FA0000, 0x37FA8000, 0x37FB0000, 0x37FB8000, 0x37FC0000, 0x37FC8000, 0x37FD0000, 0x37FD8000, 0x37FE0000, 0x37FE8000, 0x37FF0000, 0x37FF8000, 
1144				0x38000000, 0x38004000, 0x38008000, 0x3800C000, 0x38010000, 0x38014000, 0x38018000, 0x3801C000, 0x38020000, 0x38024000, 0x38028000, 0x3802C000, 0x38030000, 0x38034000, 0x38038000, 0x3803C000, 
1145				0x38040000, 0x38044000, 0x38048000, 0x3804C000, 0x38050000, 0x38054000, 0x38058000, 0x3805C000, 0x38060000, 0x38064000, 0x38068000, 0x3806C000, 0x38070000, 0x38074000, 0x38078000, 0x3807C000, 
1146				0x38080000, 0x38084000, 0x38088000, 0x3808C000, 0x38090000, 0x38094000, 0x38098000, 0x3809C000, 0x380A0000, 0x380A4000, 0x380A8000, 0x380AC000, 0x380B0000, 0x380B4000, 0x380B8000, 0x380BC000, 
1147				0x380C0000, 0x380C4000, 0x380C8000, 0x380CC000, 0x380D0000, 0x380D4000, 0x380D8000, 0x380DC000, 0x380E0000, 0x380E4000, 0x380E8000, 0x380EC000, 0x380F0000, 0x380F4000, 0x380F8000, 0x380FC000, 
1148				0x38100000, 0x38104000, 0x38108000, 0x3810C000, 0x38110000, 0x38114000, 0x38118000, 0x3811C000, 0x38120000, 0x38124000, 0x38128000, 0x3812C000, 0x38130000, 0x38134000, 0x38138000, 0x3813C000, 
1149				0x38140000, 0x38144000, 0x38148000, 0x3814C000, 0x38150000, 0x38154000, 0x38158000, 0x3815C000, 0x38160000, 0x38164000, 0x38168000, 0x3816C000, 0x38170000, 0x38174000, 0x38178000, 0x3817C000, 
1150				0x38180000, 0x38184000, 0x38188000, 0x3818C000, 0x38190000, 0x38194000, 0x38198000, 0x3819C000, 0x381A0000, 0x381A4000, 0x381A8000, 0x381AC000, 0x381B0000, 0x381B4000, 0x381B8000, 0x381BC000, 
1151				0x381C0000, 0x381C4000, 0x381C8000, 0x381CC000, 0x381D0000, 0x381D4000, 0x381D8000, 0x381DC000, 0x381E0000, 0x381E4000, 0x381E8000, 0x381EC000, 0x381F0000, 0x381F4000, 0x381F8000, 0x381FC000, 
1152				0x38200000, 0x38204000, 0x38208000, 0x3820C000, 0x38210000, 0x38214000, 0x38218000, 0x3821C000, 0x38220000, 0x38224000, 0x38228000, 0x3822C000, 0x38230000, 0x38234000, 0x38238000, 0x3823C000, 
1153				0x38240000, 0x38244000, 0x38248000, 0x3824C000, 0x38250000, 0x38254000, 0x38258000, 0x3825C000, 0x38260000, 0x38264000, 0x38268000, 0x3826C000, 0x38270000, 0x38274000, 0x38278000, 0x3827C000, 
1154				0x38280000, 0x38284000, 0x38288000, 0x3828C000, 0x38290000, 0x38294000, 0x38298000, 0x3829C000, 0x382A0000, 0x382A4000, 0x382A8000, 0x382AC000, 0x382B0000, 0x382B4000, 0x382B8000, 0x382BC000, 
1155				0x382C0000, 0x382C4000, 0x382C8000, 0x382CC000, 0x382D0000, 0x382D4000, 0x382D8000, 0x382DC000, 0x382E0000, 0x382E4000, 0x382E8000, 0x382EC000, 0x382F0000, 0x382F4000, 0x382F8000, 0x382FC000, 
1156				0x38300000, 0x38304000, 0x38308000, 0x3830C000, 0x38310000, 0x38314000, 0x38318000, 0x3831C000, 0x38320000, 0x38324000, 0x38328000, 0x3832C000, 0x38330000, 0x38334000, 0x38338000, 0x3833C000, 
1157				0x38340000, 0x38344000, 0x38348000, 0x3834C000, 0x38350000, 0x38354000, 0x38358000, 0x3835C000, 0x38360000, 0x38364000, 0x38368000, 0x3836C000, 0x38370000, 0x38374000, 0x38378000, 0x3837C000, 
1158				0x38380000, 0x38384000, 0x38388000, 0x3838C000, 0x38390000, 0x38394000, 0x38398000, 0x3839C000, 0x383A0000, 0x383A4000, 0x383A8000, 0x383AC000, 0x383B0000, 0x383B4000, 0x383B8000, 0x383BC000, 
1159				0x383C0000, 0x383C4000, 0x383C8000, 0x383CC000, 0x383D0000, 0x383D4000, 0x383D8000, 0x383DC000, 0x383E0000, 0x383E4000, 0x383E8000, 0x383EC000, 0x383F0000, 0x383F4000, 0x383F8000, 0x383FC000, 
1160				0x38400000, 0x38404000, 0x38408000, 0x3840C000, 0x38410000, 0x38414000, 0x38418000, 0x3841C000, 0x38420000, 0x38424000, 0x38428000, 0x3842C000, 0x38430000, 0x38434000, 0x38438000, 0x3843C000, 
1161				0x38440000, 0x38444000, 0x38448000, 0x3844C000, 0x38450000, 0x38454000, 0x38458000, 0x3845C000, 0x38460000, 0x38464000, 0x38468000, 0x3846C000, 0x38470000, 0x38474000, 0x38478000, 0x3847C000, 
1162				0x38480000, 0x38484000, 0x38488000, 0x3848C000, 0x38490000, 0x38494000, 0x38498000, 0x3849C000, 0x384A0000, 0x384A4000, 0x384A8000, 0x384AC000, 0x384B0000, 0x384B4000, 0x384B8000, 0x384BC000, 
1163				0x384C0000, 0x384C4000, 0x384C8000, 0x384CC000, 0x384D0000, 0x384D4000, 0x384D8000, 0x384DC000, 0x384E0000, 0x384E4000, 0x384E8000, 0x384EC000, 0x384F0000, 0x384F4000, 0x384F8000, 0x384FC000, 
1164				0x38500000, 0x38504000, 0x38508000, 0x3850C000, 0x38510000, 0x38514000, 0x38518000, 0x3851C000, 0x38520000, 0x38524000, 0x38528000, 0x3852C000, 0x38530000, 0x38534000, 0x38538000, 0x3853C000, 
1165				0x38540000, 0x38544000, 0x38548000, 0x3854C000, 0x38550000, 0x38554000, 0x38558000, 0x3855C000, 0x38560000, 0x38564000, 0x38568000, 0x3856C000, 0x38570000, 0x38574000, 0x38578000, 0x3857C000, 
1166				0x38580000, 0x38584000, 0x38588000, 0x3858C000, 0x38590000, 0x38594000, 0x38598000, 0x3859C000, 0x385A0000, 0x385A4000, 0x385A8000, 0x385AC000, 0x385B0000, 0x385B4000, 0x385B8000, 0x385BC000, 
1167				0x385C0000, 0x385C4000, 0x385C8000, 0x385CC000, 0x385D0000, 0x385D4000, 0x385D8000, 0x385DC000, 0x385E0000, 0x385E4000, 0x385E8000, 0x385EC000, 0x385F0000, 0x385F4000, 0x385F8000, 0x385FC000, 
1168				0x38600000, 0x38604000, 0x38608000, 0x3860C000, 0x38610000, 0x38614000, 0x38618000, 0x3861C000, 0x38620000, 0x38624000, 0x38628000, 0x3862C000, 0x38630000, 0x38634000, 0x38638000, 0x3863C000, 
1169				0x38640000, 0x38644000, 0x38648000, 0x3864C000, 0x38650000, 0x38654000, 0x38658000, 0x3865C000, 0x38660000, 0x38664000, 0x38668000, 0x3866C000, 0x38670000, 0x38674000, 0x38678000, 0x3867C000, 
1170				0x38680000, 0x38684000, 0x38688000, 0x3868C000, 0x38690000, 0x38694000, 0x38698000, 0x3869C000, 0x386A0000, 0x386A4000, 0x386A8000, 0x386AC000, 0x386B0000, 0x386B4000, 0x386B8000, 0x386BC000, 
1171				0x386C0000, 0x386C4000, 0x386C8000, 0x386CC000, 0x386D0000, 0x386D4000, 0x386D8000, 0x386DC000, 0x386E0000, 0x386E4000, 0x386E8000, 0x386EC000, 0x386F0000, 0x386F4000, 0x386F8000, 0x386FC000, 
1172				0x38700000, 0x38704000, 0x38708000, 0x3870C000, 0x38710000, 0x38714000, 0x38718000, 0x3871C000, 0x38720000, 0x38724000, 0x38728000, 0x3872C000, 0x38730000, 0x38734000, 0x38738000, 0x3873C000, 
1173				0x38740000, 0x38744000, 0x38748000, 0x3874C000, 0x38750000, 0x38754000, 0x38758000, 0x3875C000, 0x38760000, 0x38764000, 0x38768000, 0x3876C000, 0x38770000, 0x38774000, 0x38778000, 0x3877C000, 
1174				0x38780000, 0x38784000, 0x38788000, 0x3878C000, 0x38790000, 0x38794000, 0x38798000, 0x3879C000, 0x387A0000, 0x387A4000, 0x387A8000, 0x387AC000, 0x387B0000, 0x387B4000, 0x387B8000, 0x387BC000, 
1175				0x387C0000, 0x387C4000, 0x387C8000, 0x387CC000, 0x387D0000, 0x387D4000, 0x387D8000, 0x387DC000, 0x387E0000, 0x387E4000, 0x387E8000, 0x387EC000, 0x387F0000, 0x387F4000, 0x387F8000, 0x387FC000, 
1176				0x38000000, 0x38002000, 0x38004000, 0x38006000, 0x38008000, 0x3800A000, 0x3800C000, 0x3800E000, 0x38010000, 0x38012000, 0x38014000, 0x38016000, 0x38018000, 0x3801A000, 0x3801C000, 0x3801E000, 
1177				0x38020000, 0x38022000, 0x38024000, 0x38026000, 0x38028000, 0x3802A000, 0x3802C000, 0x3802E000, 0x38030000, 0x38032000, 0x38034000, 0x38036000, 0x38038000, 0x3803A000, 0x3803C000, 0x3803E000, 
1178				0x38040000, 0x38042000, 0x38044000, 0x38046000, 0x38048000, 0x3804A000, 0x3804C000, 0x3804E000, 0x38050000, 0x38052000, 0x38054000, 0x38056000, 0x38058000, 0x3805A000, 0x3805C000, 0x3805E000, 
1179				0x38060000, 0x38062000, 0x38064000, 0x38066000, 0x38068000, 0x3806A000, 0x3806C000, 0x3806E000, 0x38070000, 0x38072000, 0x38074000, 0x38076000, 0x38078000, 0x3807A000, 0x3807C000, 0x3807E000, 
1180				0x38080000, 0x38082000, 0x38084000, 0x38086000, 0x38088000, 0x3808A000, 0x3808C000, 0x3808E000, 0x38090000, 0x38092000, 0x38094000, 0x38096000, 0x38098000, 0x3809A000, 0x3809C000, 0x3809E000, 
1181				0x380A0000, 0x380A2000, 0x380A4000, 0x380A6000, 0x380A8000, 0x380AA000, 0x380AC000, 0x380AE000, 0x380B0000, 0x380B2000, 0x380B4000, 0x380B6000, 0x380B8000, 0x380BA000, 0x380BC000, 0x380BE000, 
1182				0x380C0000, 0x380C2000, 0x380C4000, 0x380C6000, 0x380C8000, 0x380CA000, 0x380CC000, 0x380CE000, 0x380D0000, 0x380D2000, 0x380D4000, 0x380D6000, 0x380D8000, 0x380DA000, 0x380DC000, 0x380DE000, 
1183				0x380E0000, 0x380E2000, 0x380E4000, 0x380E6000, 0x380E8000, 0x380EA000, 0x380EC000, 0x380EE000, 0x380F0000, 0x380F2000, 0x380F4000, 0x380F6000, 0x380F8000, 0x380FA000, 0x380FC000, 0x380FE000, 
1184				0x38100000, 0x38102000, 0x38104000, 0x38106000, 0x38108000, 0x3810A000, 0x3810C000, 0x3810E000, 0x38110000, 0x38112000, 0x38114000, 0x38116000, 0x38118000, 0x3811A000, 0x3811C000, 0x3811E000, 
1185				0x38120000, 0x38122000, 0x38124000, 0x38126000, 0x38128000, 0x3812A000, 0x3812C000, 0x3812E000, 0x38130000, 0x38132000, 0x38134000, 0x38136000, 0x38138000, 0x3813A000, 0x3813C000, 0x3813E000, 
1186				0x38140000, 0x38142000, 0x38144000, 0x38146000, 0x38148000, 0x3814A000, 0x3814C000, 0x3814E000, 0x38150000, 0x38152000, 0x38154000, 0x38156000, 0x38158000, 0x3815A000, 0x3815C000, 0x3815E000, 
1187				0x38160000, 0x38162000, 0x38164000, 0x38166000, 0x38168000, 0x3816A000, 0x3816C000, 0x3816E000, 0x38170000, 0x38172000, 0x38174000, 0x38176000, 0x38178000, 0x3817A000, 0x3817C000, 0x3817E000, 
1188				0x38180000, 0x38182000, 0x38184000, 0x38186000, 0x38188000, 0x3818A000, 0x3818C000, 0x3818E000, 0x38190000, 0x38192000, 0x38194000, 0x38196000, 0x38198000, 0x3819A000, 0x3819C000, 0x3819E000, 
1189				0x381A0000, 0x381A2000, 0x381A4000, 0x381A6000, 0x381A8000, 0x381AA000, 0x381AC000, 0x381AE000, 0x381B0000, 0x381B2000, 0x381B4000, 0x381B6000, 0x381B8000, 0x381BA000, 0x381BC000, 0x381BE000, 
1190				0x381C0000, 0x381C2000, 0x381C4000, 0x381C6000, 0x381C8000, 0x381CA000, 0x381CC000, 0x381CE000, 0x381D0000, 0x381D2000, 0x381D4000, 0x381D6000, 0x381D8000, 0x381DA000, 0x381DC000, 0x381DE000, 
1191				0x381E0000, 0x381E2000, 0x381E4000, 0x381E6000, 0x381E8000, 0x381EA000, 0x381EC000, 0x381EE000, 0x381F0000, 0x381F2000, 0x381F4000, 0x381F6000, 0x381F8000, 0x381FA000, 0x381FC000, 0x381FE000, 
1192				0x38200000, 0x38202000, 0x38204000, 0x38206000, 0x38208000, 0x3820A000, 0x3820C000, 0x3820E000, 0x38210000, 0x38212000, 0x38214000, 0x38216000, 0x38218000, 0x3821A000, 0x3821C000, 0x3821E000, 
1193				0x38220000, 0x38222000, 0x38224000, 0x38226000, 0x38228000, 0x3822A000, 0x3822C000, 0x3822E000, 0x38230000, 0x38232000, 0x38234000, 0x38236000, 0x38238000, 0x3823A000, 0x3823C000, 0x3823E000, 
1194				0x38240000, 0x38242000, 0x38244000, 0x38246000, 0x38248000, 0x3824A000, 0x3824C000, 0x3824E000, 0x38250000, 0x38252000, 0x38254000, 0x38256000, 0x38258000, 0x3825A000, 0x3825C000, 0x3825E000, 
1195				0x38260000, 0x38262000, 0x38264000, 0x38266000, 0x38268000, 0x3826A000, 0x3826C000, 0x3826E000, 0x38270000, 0x38272000, 0x38274000, 0x38276000, 0x38278000, 0x3827A000, 0x3827C000, 0x3827E000, 
1196				0x38280000, 0x38282000, 0x38284000, 0x38286000, 0x38288000, 0x3828A000, 0x3828C000, 0x3828E000, 0x38290000, 0x38292000, 0x38294000, 0x38296000, 0x38298000, 0x3829A000, 0x3829C000, 0x3829E000, 
1197				0x382A0000, 0x382A2000, 0x382A4000, 0x382A6000, 0x382A8000, 0x382AA000, 0x382AC000, 0x382AE000, 0x382B0000, 0x382B2000, 0x382B4000, 0x382B6000, 0x382B8000, 0x382BA000, 0x382BC000, 0x382BE000, 
1198				0x382C0000, 0x382C2000, 0x382C4000, 0x382C6000, 0x382C8000, 0x382CA000, 0x382CC000, 0x382CE000, 0x382D0000, 0x382D2000, 0x382D4000, 0x382D6000, 0x382D8000, 0x382DA000, 0x382DC000, 0x382DE000, 
1199				0x382E0000, 0x382E2000, 0x382E4000, 0x382E6000, 0x382E8000, 0x382EA000, 0x382EC000, 0x382EE000, 0x382F0000, 0x382F2000, 0x382F4000, 0x382F6000, 0x382F8000, 0x382FA000, 0x382FC000, 0x382FE000, 
1200				0x38300000, 0x38302000, 0x38304000, 0x38306000, 0x38308000, 0x3830A000, 0x3830C000, 0x3830E000, 0x38310000, 0x38312000, 0x38314000, 0x38316000, 0x38318000, 0x3831A000, 0x3831C000, 0x3831E000, 
1201				0x38320000, 0x38322000, 0x38324000, 0x38326000, 0x38328000, 0x3832A000, 0x3832C000, 0x3832E000, 0x38330000, 0x38332000, 0x38334000, 0x38336000, 0x38338000, 0x3833A000, 0x3833C000, 0x3833E000, 
1202				0x38340000, 0x38342000, 0x38344000, 0x38346000, 0x38348000, 0x3834A000, 0x3834C000, 0x3834E000, 0x38350000, 0x38352000, 0x38354000, 0x38356000, 0x38358000, 0x3835A000, 0x3835C000, 0x3835E000, 
1203				0x38360000, 0x38362000, 0x38364000, 0x38366000, 0x38368000, 0x3836A000, 0x3836C000, 0x3836E000, 0x38370000, 0x38372000, 0x38374000, 0x38376000, 0x38378000, 0x3837A000, 0x3837C000, 0x3837E000, 
1204				0x38380000, 0x38382000, 0x38384000, 0x38386000, 0x38388000, 0x3838A000, 0x3838C000, 0x3838E000, 0x38390000, 0x38392000, 0x38394000, 0x38396000, 0x38398000, 0x3839A000, 0x3839C000, 0x3839E000, 
1205				0x383A0000, 0x383A2000, 0x383A4000, 0x383A6000, 0x383A8000, 0x383AA000, 0x383AC000, 0x383AE000, 0x383B0000, 0x383B2000, 0x383B4000, 0x383B6000, 0x383B8000, 0x383BA000, 0x383BC000, 0x383BE000, 
1206				0x383C0000, 0x383C2000, 0x383C4000, 0x383C6000, 0x383C8000, 0x383CA000, 0x383CC000, 0x383CE000, 0x383D0000, 0x383D2000, 0x383D4000, 0x383D6000, 0x383D8000, 0x383DA000, 0x383DC000, 0x383DE000, 
1207				0x383E0000, 0x383E2000, 0x383E4000, 0x383E6000, 0x383E8000, 0x383EA000, 0x383EC000, 0x383EE000, 0x383F0000, 0x383F2000, 0x383F4000, 0x383F6000, 0x383F8000, 0x383FA000, 0x383FC000, 0x383FE000, 
1208				0x38400000, 0x38402000, 0x38404000, 0x38406000, 0x38408000, 0x3840A000, 0x3840C000, 0x3840E000, 0x38410000, 0x38412000, 0x38414000, 0x38416000, 0x38418000, 0x3841A000, 0x3841C000, 0x3841E000, 
1209				0x38420000, 0x38422000, 0x38424000, 0x38426000, 0x38428000, 0x3842A000, 0x3842C000, 0x3842E000, 0x38430000, 0x38432000, 0x38434000, 0x38436000, 0x38438000, 0x3843A000, 0x3843C000, 0x3843E000, 
1210				0x38440000, 0x38442000, 0x38444000, 0x38446000, 0x38448000, 0x3844A000, 0x3844C000, 0x3844E000, 0x38450000, 0x38452000, 0x38454000, 0x38456000, 0x38458000, 0x3845A000, 0x3845C000, 0x3845E000, 
1211				0x38460000, 0x38462000, 0x38464000, 0x38466000, 0x38468000, 0x3846A000, 0x3846C000, 0x3846E000, 0x38470000, 0x38472000, 0x38474000, 0x38476000, 0x38478000, 0x3847A000, 0x3847C000, 0x3847E000, 
1212				0x38480000, 0x38482000, 0x38484000, 0x38486000, 0x38488000, 0x3848A000, 0x3848C000, 0x3848E000, 0x38490000, 0x38492000, 0x38494000, 0x38496000, 0x38498000, 0x3849A000, 0x3849C000, 0x3849E000, 
1213				0x384A0000, 0x384A2000, 0x384A4000, 0x384A6000, 0x384A8000, 0x384AA000, 0x384AC000, 0x384AE000, 0x384B0000, 0x384B2000, 0x384B4000, 0x384B6000, 0x384B8000, 0x384BA000, 0x384BC000, 0x384BE000, 
1214				0x384C0000, 0x384C2000, 0x384C4000, 0x384C6000, 0x384C8000, 0x384CA000, 0x384CC000, 0x384CE000, 0x384D0000, 0x384D2000, 0x384D4000, 0x384D6000, 0x384D8000, 0x384DA000, 0x384DC000, 0x384DE000, 
1215				0x384E0000, 0x384E2000, 0x384E4000, 0x384E6000, 0x384E8000, 0x384EA000, 0x384EC000, 0x384EE000, 0x384F0000, 0x384F2000, 0x384F4000, 0x384F6000, 0x384F8000, 0x384FA000, 0x384FC000, 0x384FE000, 
1216				0x38500000, 0x38502000, 0x38504000, 0x38506000, 0x38508000, 0x3850A000, 0x3850C000, 0x3850E000, 0x38510000, 0x38512000, 0x38514000, 0x38516000, 0x38518000, 0x3851A000, 0x3851C000, 0x3851E000, 
1217				0x38520000, 0x38522000, 0x38524000, 0x38526000, 0x38528000, 0x3852A000, 0x3852C000, 0x3852E000, 0x38530000, 0x38532000, 0x38534000, 0x38536000, 0x38538000, 0x3853A000, 0x3853C000, 0x3853E000, 
1218				0x38540000, 0x38542000, 0x38544000, 0x38546000, 0x38548000, 0x3854A000, 0x3854C000, 0x3854E000, 0x38550000, 0x38552000, 0x38554000, 0x38556000, 0x38558000, 0x3855A000, 0x3855C000, 0x3855E000, 
1219				0x38560000, 0x38562000, 0x38564000, 0x38566000, 0x38568000, 0x3856A000, 0x3856C000, 0x3856E000, 0x38570000, 0x38572000, 0x38574000, 0x38576000, 0x38578000, 0x3857A000, 0x3857C000, 0x3857E000, 
1220				0x38580000, 0x38582000, 0x38584000, 0x38586000, 0x38588000, 0x3858A000, 0x3858C000, 0x3858E000, 0x38590000, 0x38592000, 0x38594000, 0x38596000, 0x38598000, 0x3859A000, 0x3859C000, 0x3859E000, 
1221				0x385A0000, 0x385A2000, 0x385A4000, 0x385A6000, 0x385A8000, 0x385AA000, 0x385AC000, 0x385AE000, 0x385B0000, 0x385B2000, 0x385B4000, 0x385B6000, 0x385B8000, 0x385BA000, 0x385BC000, 0x385BE000, 
1222				0x385C0000, 0x385C2000, 0x385C4000, 0x385C6000, 0x385C8000, 0x385CA000, 0x385CC000, 0x385CE000, 0x385D0000, 0x385D2000, 0x385D4000, 0x385D6000, 0x385D8000, 0x385DA000, 0x385DC000, 0x385DE000, 
1223				0x385E0000, 0x385E2000, 0x385E4000, 0x385E6000, 0x385E8000, 0x385EA000, 0x385EC000, 0x385EE000, 0x385F0000, 0x385F2000, 0x385F4000, 0x385F6000, 0x385F8000, 0x385FA000, 0x385FC000, 0x385FE000, 
1224				0x38600000, 0x38602000, 0x38604000, 0x38606000, 0x38608000, 0x3860A000, 0x3860C000, 0x3860E000, 0x38610000, 0x38612000, 0x38614000, 0x38616000, 0x38618000, 0x3861A000, 0x3861C000, 0x3861E000, 
1225				0x38620000, 0x38622000, 0x38624000, 0x38626000, 0x38628000, 0x3862A000, 0x3862C000, 0x3862E000, 0x38630000, 0x38632000, 0x38634000, 0x38636000, 0x38638000, 0x3863A000, 0x3863C000, 0x3863E000, 
1226				0x38640000, 0x38642000, 0x38644000, 0x38646000, 0x38648000, 0x3864A000, 0x3864C000, 0x3864E000, 0x38650000, 0x38652000, 0x38654000, 0x38656000, 0x38658000, 0x3865A000, 0x3865C000, 0x3865E000, 
1227				0x38660000, 0x38662000, 0x38664000, 0x38666000, 0x38668000, 0x3866A000, 0x3866C000, 0x3866E000, 0x38670000, 0x38672000, 0x38674000, 0x38676000, 0x38678000, 0x3867A000, 0x3867C000, 0x3867E000, 
1228				0x38680000, 0x38682000, 0x38684000, 0x38686000, 0x38688000, 0x3868A000, 0x3868C000, 0x3868E000, 0x38690000, 0x38692000, 0x38694000, 0x38696000, 0x38698000, 0x3869A000, 0x3869C000, 0x3869E000, 
1229				0x386A0000, 0x386A2000, 0x386A4000, 0x386A6000, 0x386A8000, 0x386AA000, 0x386AC000, 0x386AE000, 0x386B0000, 0x386B2000, 0x386B4000, 0x386B6000, 0x386B8000, 0x386BA000, 0x386BC000, 0x386BE000, 
1230				0x386C0000, 0x386C2000, 0x386C4000, 0x386C6000, 0x386C8000, 0x386CA000, 0x386CC000, 0x386CE000, 0x386D0000, 0x386D2000, 0x386D4000, 0x386D6000, 0x386D8000, 0x386DA000, 0x386DC000, 0x386DE000, 
1231				0x386E0000, 0x386E2000, 0x386E4000, 0x386E6000, 0x386E8000, 0x386EA000, 0x386EC000, 0x386EE000, 0x386F0000, 0x386F2000, 0x386F4000, 0x386F6000, 0x386F8000, 0x386FA000, 0x386FC000, 0x386FE000, 
1232				0x38700000, 0x38702000, 0x38704000, 0x38706000, 0x38708000, 0x3870A000, 0x3870C000, 0x3870E000, 0x38710000, 0x38712000, 0x38714000, 0x38716000, 0x38718000, 0x3871A000, 0x3871C000, 0x3871E000, 
1233				0x38720000, 0x38722000, 0x38724000, 0x38726000, 0x38728000, 0x3872A000, 0x3872C000, 0x3872E000, 0x38730000, 0x38732000, 0x38734000, 0x38736000, 0x38738000, 0x3873A000, 0x3873C000, 0x3873E000, 
1234				0x38740000, 0x38742000, 0x38744000, 0x38746000, 0x38748000, 0x3874A000, 0x3874C000, 0x3874E000, 0x38750000, 0x38752000, 0x38754000, 0x38756000, 0x38758000, 0x3875A000, 0x3875C000, 0x3875E000, 
1235				0x38760000, 0x38762000, 0x38764000, 0x38766000, 0x38768000, 0x3876A000, 0x3876C000, 0x3876E000, 0x38770000, 0x38772000, 0x38774000, 0x38776000, 0x38778000, 0x3877A000, 0x3877C000, 0x3877E000, 
1236				0x38780000, 0x38782000, 0x38784000, 0x38786000, 0x38788000, 0x3878A000, 0x3878C000, 0x3878E000, 0x38790000, 0x38792000, 0x38794000, 0x38796000, 0x38798000, 0x3879A000, 0x3879C000, 0x3879E000, 
1237				0x387A0000, 0x387A2000, 0x387A4000, 0x387A6000, 0x387A8000, 0x387AA000, 0x387AC000, 0x387AE000, 0x387B0000, 0x387B2000, 0x387B4000, 0x387B6000, 0x387B8000, 0x387BA000, 0x387BC000, 0x387BE000, 
1238				0x387C0000, 0x387C2000, 0x387C4000, 0x387C6000, 0x387C8000, 0x387CA000, 0x387CC000, 0x387CE000, 0x387D0000, 0x387D2000, 0x387D4000, 0x387D6000, 0x387D8000, 0x387DA000, 0x387DC000, 0x387DE000, 
1239				0x387E0000, 0x387E2000, 0x387E4000, 0x387E6000, 0x387E8000, 0x387EA000, 0x387EC000, 0x387EE000, 0x387F0000, 0x387F2000, 0x387F4000, 0x387F6000, 0x387F8000, 0x387FA000, 0x387FC000, 0x387FE000 };
1240			static const bits<float>::type exponent_table[64] = {
1241				0x00000000, 0x00800000, 0x01000000, 0x01800000, 0x02000000, 0x02800000, 0x03000000, 0x03800000, 0x04000000, 0x04800000, 0x05000000, 0x05800000, 0x06000000, 0x06800000, 0x07000000, 0x07800000, 
1242				0x08000000, 0x08800000, 0x09000000, 0x09800000, 0x0A000000, 0x0A800000, 0x0B000000, 0x0B800000, 0x0C000000, 0x0C800000, 0x0D000000, 0x0D800000, 0x0E000000, 0x0E800000, 0x0F000000, 0x47800000, 
1243				0x80000000, 0x80800000, 0x81000000, 0x81800000, 0x82000000, 0x82800000, 0x83000000, 0x83800000, 0x84000000, 0x84800000, 0x85000000, 0x85800000, 0x86000000, 0x86800000, 0x87000000, 0x87800000, 
1244				0x88000000, 0x88800000, 0x89000000, 0x89800000, 0x8A000000, 0x8A800000, 0x8B000000, 0x8B800000, 0x8C000000, 0x8C800000, 0x8D000000, 0x8D800000, 0x8E000000, 0x8E800000, 0x8F000000, 0xC7800000 };
1245			static const unsigned short offset_table[64] = {
1246				0, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 
1247				0, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024, 1024 };
1248			bits<float>::type fbits = mantissa_table[offset_table[value>>10]+(value&0x3FF)] + exponent_table[value>>10];
1249		#endif
1250			float out;
1251			std::memcpy(&out, &fbits, sizeof(float));
1252			return out;
1253		#endif
1254		}
1255
1256		/// Convert half-precision to IEEE double-precision.
1257		/// \param value half-precision value to convert
1258		/// \return double-precision value
1259		inline double half2float_impl(unsigned int value, double, true_type)
1260		{
1261		#if HALF_ENABLE_F16C_INTRINSICS
1262			return _mm_cvtsd_f64(_mm_cvtps_pd(_mm_cvtph_ps(_mm_cvtsi32_si128(value))));
1263		#else
1264			uint32 hi = static_cast<uint32>(value&0x8000) << 16;
1265			unsigned int abs = value & 0x7FFF;
1266			if(abs)
1267			{
1268				hi |= 0x3F000000 << static_cast<unsigned>(abs>=0x7C00);
1269				for(; abs<0x400; abs<<=1,hi-=0x100000) ;
1270				hi += static_cast<uint32>(abs) << 10;
1271			}
1272			bits<double>::type dbits = static_cast<bits<double>::type>(hi) << 32;
1273			double out;
1274			std::memcpy(&out, &dbits, sizeof(double));
1275			return out;
1276		#endif
1277		}
1278
1279		/// Convert half-precision to non-IEEE floating-point.
1280		/// \tparam T type to convert to (builtin integer type)
1281		/// \param value half-precision value to convert
1282		/// \return floating-point value
1283		template<typename T> T half2float_impl(unsigned int value, T, ...)
1284		{
1285			T out;
1286			unsigned int abs = value & 0x7FFF;
1287			if(abs > 0x7C00)
1288				out = (std::numeric_limits<T>::has_signaling_NaN && !(abs&0x200)) ? std::numeric_limits<T>::signaling_NaN() :
1289					std::numeric_limits<T>::has_quiet_NaN ? std::numeric_limits<T>::quiet_NaN() : T();
1290			else if(abs == 0x7C00)
1291				out = std::numeric_limits<T>::has_infinity ? std::numeric_limits<T>::infinity() : std::numeric_limits<T>::max();
1292			else if(abs > 0x3FF)
1293				out = std::ldexp(static_cast<T>((abs&0x3FF)|0x400), (abs>>10)-25);
1294			else
1295				out = std::ldexp(static_cast<T>(abs), -24);
1296			return (value&0x8000) ? -out : out;
1297		}
1298
1299		/// Convert half-precision to floating-point.
1300		/// \tparam T type to convert to (builtin integer type)
1301		/// \param value half-precision value to convert
1302		/// \return floating-point value
1303		template<typename T> T half2float(unsigned int value)
1304		{
1305			return half2float_impl(value, T(), bool_type<std::numeric_limits<T>::is_iec559&&sizeof(typename bits<T>::type)==sizeof(T)>());
1306		}
1307
1308		/// Convert half-precision floating-point to integer.
1309		/// \tparam R rounding mode to use
1310		/// \tparam E `true` for round to even, `false` for round away from zero
1311		/// \tparam I `true` to raise INEXACT exception (if inexact), `false` to never raise it
1312		/// \tparam T type to convert to (buitlin integer type with at least 16 bits precision, excluding any implicit sign bits)
1313		/// \param value half-precision value to convert
1314		/// \return rounded integer value
1315		/// \exception FE_INVALID if value is not representable in type \a T
1316		/// \exception FE_INEXACT if value had to be rounded and \a I is `true`
1317		template<std::float_round_style R,bool E,bool I,typename T> T half2int(unsigned int value)
1318		{
1319			unsigned int abs = value & 0x7FFF;
1320			if(abs >= 0x7C00)
1321			{
1322				raise(FE_INVALID);
1323				return (value&0x8000) ? std::numeric_limits<T>::min() : std::numeric_limits<T>::max();
1324			}
1325			if(abs < 0x3800)
1326			{
1327				raise(FE_INEXACT, I);
1328				return	(R==std::round_toward_infinity) ? T(~(value>>15)&(abs!=0)) :
1329						(R==std::round_toward_neg_infinity) ? -T(value>0x8000) :
1330						T();
1331			}
1332			int exp = 25 - (abs>>10);
1333			unsigned int m = (value&0x3FF) | 0x400;
1334			int32 i = static_cast<int32>((exp<=0) ? (m<<-exp) : ((m+(
1335				(R==std::round_to_nearest) ? ((1<<(exp-1))-(~(m>>exp)&E)) :
1336				(R==std::round_toward_infinity) ? (((1<<exp)-1)&((value>>15)-1)) :
1337				(R==std::round_toward_neg_infinity) ? (((1<<exp)-1)&-(value>>15)) : 0))>>exp));
1338			if((!std::numeric_limits<T>::is_signed && (value&0x8000)) || (std::numeric_limits<T>::digits<16 &&
1339				((value&0x8000) ? (-i<std::numeric_limits<T>::min()) : (i>std::numeric_limits<T>::max()))))
1340				raise(FE_INVALID);
1341			else if(I && exp > 0 && (m&((1<<exp)-1)))
1342				raise(FE_INEXACT);
1343			return static_cast<T>((value&0x8000) ? -i : i);
1344		}
1345
1346		/// \}
1347		/// \name Mathematics
1348		/// \{
1349
1350		/// upper part of 64-bit multiplication.
1351		/// \tparam R rounding mode to use
1352		/// \param x first factor
1353		/// \param y second factor
1354		/// \return upper 32 bit of \a x * \a y
1355		template<std::float_round_style R> uint32 mulhi(uint32 x, uint32 y)
1356		{
1357			uint32 xy = (x>>16) * (y&0xFFFF), yx = (x&0xFFFF) * (y>>16), c = (xy&0xFFFF) + (yx&0xFFFF) + (((x&0xFFFF)*(y&0xFFFF))>>16);
1358			return (x>>16)*(y>>16) + (xy>>16) + (yx>>16) + (c>>16) +
1359				((R==std::round_to_nearest) ? ((c>>15)&1) : (R==std::round_toward_infinity) ? ((c&0xFFFF)!=0) : 0);
1360		}
1361
1362		/// 64-bit multiplication.
1363		/// \param x first factor
1364		/// \param y second factor
1365		/// \return upper 32 bit of \a x * \a y rounded to nearest
1366		inline uint32 multiply64(uint32 x, uint32 y)
1367		{
1368		#if HALF_ENABLE_CPP11_LONG_LONG
1369			return static_cast<uint32>((static_cast<unsigned long long>(x)*static_cast<unsigned long long>(y)+0x80000000)>>32);
1370		#else
1371			return mulhi<std::round_to_nearest>(x, y);
1372		#endif
1373		}
1374
1375		/// 64-bit division.
1376		/// \param x upper 32 bit of dividend
1377		/// \param y divisor
1378		/// \param s variable to store sticky bit for rounding
1379		/// \return (\a x << 32) / \a y
1380		inline uint32 divide64(uint32 x, uint32 y, int &s)
1381		{
1382		#if HALF_ENABLE_CPP11_LONG_LONG
1383			unsigned long long xx = static_cast<unsigned long long>(x) << 32;
1384			return s = (xx%y!=0), static_cast<uint32>(xx/y);
1385		#else
1386			y >>= 1;
1387			uint32 rem = x, div = 0;
1388			for(unsigned int i=0; i<32; ++i)
1389			{
1390				div <<= 1;
1391				if(rem >= y)
1392				{
1393					rem -= y;
1394					div |= 1;
1395				}
1396				rem <<= 1;
1397			}
1398			return s = rem > 1, div;
1399		#endif
1400		}
1401
1402		/// Half precision positive modulus.
1403		/// \tparam Q `true` to compute full quotient, `false` else
1404		/// \tparam R `true` to compute signed remainder, `false` for positive remainder
1405		/// \param x first operand as positive finite half-precision value
1406		/// \param y second operand as positive finite half-precision value
1407		/// \param quo adress to store quotient at, `nullptr` if \a Q `false`
1408		/// \return modulus of \a x / \a y
1409		template<bool Q,bool R> unsigned int mod(unsigned int x, unsigned int y, int *quo = NULL)
1410		{
1411			unsigned int q = 0;
1412			if(x > y)
1413			{
1414				int absx = x, absy = y, expx = 0, expy = 0;
1415				for(; absx<0x400; absx<<=1,--expx) ;
1416				for(; absy<0x400; absy<<=1,--expy) ;
1417				expx += absx >> 10;
1418				expy += absy >> 10;
1419				int mx = (absx&0x3FF) | 0x400, my = (absy&0x3FF) | 0x400;
1420				for(int d=expx-expy; d; --d)
1421				{
1422					if(!Q && mx == my)
1423						return 0;
1424					if(mx >= my)
1425					{
1426						mx -= my;
1427						q += Q;
1428					}
1429					mx <<= 1;
1430					q <<= static_cast<int>(Q);
1431				}
1432				if(!Q && mx == my)
1433					return 0;
1434				if(mx >= my)
1435				{
1436					mx -= my;
1437					++q;
1438				}
1439				if(Q)
1440				{
1441					q &= (1<<(std::numeric_limits<int>::digits-1)) - 1;
1442					if(!mx)
1443						return *quo = q, 0;
1444				}
1445				for(; mx<0x400; mx<<=1,--expy) ;
1446				x = (expy>0) ? ((expy<<10)|(mx&0x3FF)) : (mx>>(1-expy));
1447			}
1448			if(R)
1449			{
1450				unsigned int a, b;
1451				if(y < 0x800)
1452				{
1453					a = (x<0x400) ? (x<<1) : (x+0x400);
1454					b = y;
1455				}
1456				else
1457				{
1458					a = x;
1459					b = y - 0x400;
1460				}
1461				if(a > b || (a == b && (q&1)))
1462				{
1463					int exp = (y>>10) + (y<=0x3FF), d = exp - (x>>10) - (x<=0x3FF);
1464					int m = (((y&0x3FF)|((y>0x3FF)<<10))<<1) - (((x&0x3FF)|((x>0x3FF)<<10))<<(1-d));
1465					for(; m<0x800 && exp>1; m<<=1,--exp) ;
1466					x = 0x8000 + ((exp-1)<<10) + (m>>1);
1467					q += Q;
1468				}
1469			}
1470			if(Q)
1471				*quo = q;
1472			return x;
1473		}
1474
1475		/// Fixed point square root.
1476		/// \tparam F number of fractional bits
1477		/// \param r radicand in Q1.F fixed point format
1478		/// \param exp exponent
1479		/// \return square root as Q1.F/2
1480		template<unsigned int F> uint32 sqrt(uint32 &r, int &exp)
1481		{
1482			int i = exp & 1;
1483			r <<= i;
1484			exp = (exp-i) / 2;
1485			uint32 m = 0;
1486			for(uint32 bit=static_cast<uint32>(1)<<F; bit; bit>>=2)
1487			{
1488				if(r < m+bit)
1489					m >>= 1;
1490				else
1491				{
1492					r -= m + bit;
1493					m = (m>>1) + bit;
1494				}
1495			}
1496			return m;
1497		}
1498
1499		/// Fixed point binary exponential.
1500		/// This uses the BKM algorithm in E-mode.
1501		/// \param m exponent in [0,1) as Q0.31
1502		/// \param n number of iterations (at most 32)
1503		/// \return 2 ^ \a m as Q1.31
1504		inline uint32 exp2(uint32 m, unsigned int n = 32)
1505		{
1506			static const uint32 logs[] = {
1507				0x80000000, 0x4AE00D1D, 0x2934F098, 0x15C01A3A, 0x0B31FB7D, 0x05AEB4DD, 0x02DCF2D1, 0x016FE50B,
1508				0x00B84E23, 0x005C3E10, 0x002E24CA, 0x001713D6, 0x000B8A47, 0x0005C53B, 0x0002E2A3, 0x00017153,
1509				0x0000B8AA, 0x00005C55, 0x00002E2B, 0x00001715, 0x00000B8B, 0x000005C5, 0x000002E3, 0x00000171,
1510				0x000000B9, 0x0000005C, 0x0000002E, 0x00000017, 0x0000000C, 0x00000006, 0x00000003, 0x00000001 };
1511			if(!m)
1512				return 0x80000000;
1513			uint32 mx = 0x80000000, my = 0;
1514			for(unsigned int i=1; i<n; ++i)
1515			{
1516				uint32 mz = my + logs[i];
1517				if(mz <= m)
1518				{
1519					my = mz;
1520					mx += mx >> i;
1521				}
1522			}
1523			return mx;
1524		}
1525
1526		/// Fixed point binary logarithm.
1527		/// This uses the BKM algorithm in L-mode.
1528		/// \param m mantissa in [1,2) as Q1.30
1529		/// \param n number of iterations (at most 32)
1530		/// \return log2(\a m) as Q0.31
1531		inline uint32 log2(uint32 m, unsigned int n = 32)
1532		{
1533			static const uint32 logs[] = {
1534				0x80000000, 0x4AE00D1D, 0x2934F098, 0x15C01A3A, 0x0B31FB7D, 0x05AEB4DD, 0x02DCF2D1, 0x016FE50B,
1535				0x00B84E23, 0x005C3E10, 0x002E24CA, 0x001713D6, 0x000B8A47, 0x0005C53B, 0x0002E2A3, 0x00017153,
1536				0x0000B8AA, 0x00005C55, 0x00002E2B, 0x00001715, 0x00000B8B, 0x000005C5, 0x000002E3, 0x00000171,
1537				0x000000B9, 0x0000005C, 0x0000002E, 0x00000017, 0x0000000C, 0x00000006, 0x00000003, 0x00000001 };
1538			if(m == 0x40000000)
1539				return 0;
1540			uint32 mx = 0x40000000, my = 0;
1541			for(unsigned int i=1; i<n; ++i)
1542			{
1543				uint32 mz = mx + (mx>>i);
1544				if(mz <= m)
1545				{
1546					mx = mz;
1547					my += logs[i];
1548				}
1549			}
1550			return my;
1551		}
1552
1553		/// Fixed point sine and cosine.
1554		/// This uses the CORDIC algorithm in rotation mode.
1555		/// \param mz angle in [-pi/2,pi/2] as Q1.30
1556		/// \param n number of iterations (at most 31)
1557		/// \return sine and cosine of \a mz as Q1.30
1558		inline std::pair<uint32,uint32> sincos(uint32 mz, unsigned int n = 31)
1559		{
1560			static const uint32 angles[] = {
1561				0x3243F6A9, 0x1DAC6705, 0x0FADBAFD, 0x07F56EA7, 0x03FEAB77, 0x01FFD55C, 0x00FFFAAB, 0x007FFF55,
1562				0x003FFFEB, 0x001FFFFD, 0x00100000, 0x00080000, 0x00040000, 0x00020000, 0x00010000, 0x00008000,
1563				0x00004000, 0x00002000, 0x00001000, 0x00000800, 0x00000400, 0x00000200, 0x00000100, 0x00000080,
1564				0x00000040, 0x00000020, 0x00000010, 0x00000008, 0x00000004, 0x00000002, 0x00000001 };
1565			uint32 mx = 0x26DD3B6A, my = 0;
1566			for(unsigned int i=0; i<n; ++i)
1567			{
1568				uint32 sign = sign_mask(mz);
1569				uint32 tx = mx - (arithmetic_shift(my, i)^sign) + sign;
1570				uint32 ty = my + (arithmetic_shift(mx, i)^sign) - sign;
1571				mx = tx; my = ty; mz -= (angles[i]^sign) - sign;
1572			}
1573			return std::make_pair(my, mx);
1574		}
1575
1576		/// Fixed point arc tangent.
1577		/// This uses the CORDIC algorithm in vectoring mode.
1578		/// \param my y coordinate as Q0.30
1579		/// \param mx x coordinate as Q0.30
1580		/// \param n number of iterations (at most 31)
1581		/// \return arc tangent of \a my / \a mx as Q1.30
1582		inline uint32 atan2(uint32 my, uint32 mx, unsigned int n = 31)
1583		{
1584			static const uint32 angles[] = {
1585				0x3243F6A9, 0x1DAC6705, 0x0FADBAFD, 0x07F56EA7, 0x03FEAB77, 0x01FFD55C, 0x00FFFAAB, 0x007FFF55,
1586				0x003FFFEB, 0x001FFFFD, 0x00100000, 0x00080000, 0x00040000, 0x00020000, 0x00010000, 0x00008000,
1587				0x00004000, 0x00002000, 0x00001000, 0x00000800, 0x00000400, 0x00000200, 0x00000100, 0x00000080,
1588				0x00000040, 0x00000020, 0x00000010, 0x00000008, 0x00000004, 0x00000002, 0x00000001 };
1589			uint32 mz = 0;
1590			for(unsigned int i=0; i<n; ++i)
1591			{
1592				uint32 sign = sign_mask(my);
1593				uint32 tx = mx + (arithmetic_shift(my, i)^sign) - sign;
1594				uint32 ty = my - (arithmetic_shift(mx, i)^sign) + sign;
1595				mx = tx; my = ty; mz += (angles[i]^sign) - sign;
1596			}
1597			return mz;
1598		}
1599
1600		/// Reduce argument for trigonometric functions.
1601		/// \param abs half-precision floating-point value
1602		/// \param k value to take quarter period
1603		/// \return \a abs reduced to [-pi/4,pi/4] as Q0.30
1604		inline uint32 angle_arg(unsigned int abs, int &k)
1605		{
1606			uint32 m = (abs&0x3FF) | ((abs>0x3FF)<<10);
1607			int exp = (abs>>10) + (abs<=0x3FF) - 15;
1608			if(abs < 0x3A48)
1609				return k = 0, m << (exp+20);
1610		#if HALF_ENABLE_CPP11_LONG_LONG
1611			unsigned long long y = m * 0xA2F9836E4E442, mask = (1ULL<<(62-exp)) - 1, yi = (y+(mask>>1)) & ~mask, f = y - yi;
1612			uint32 sign = -static_cast<uint32>(f>>63);
1613			k = static_cast<int>(yi>>(62-exp));
1614			return (multiply64(static_cast<uint32>((sign ? -f : f)>>(31-exp)), 0xC90FDAA2)^sign) - sign;
1615		#else
1616			uint32 yh = m*0xA2F98 + mulhi<std::round_toward_zero>(m, 0x36E4E442), yl = (m*0x36E4E442) & 0xFFFFFFFF;
1617			uint32 mask = (static_cast<uint32>(1)<<(30-exp)) - 1, yi = (yh+(mask>>1)) & ~mask, sign = -static_cast<uint32>(yi>yh);
1618			k = static_cast<int>(yi>>(30-exp));
1619			uint32 fh = (yh^sign) + (yi^~sign) - ~sign, fl = (yl^sign) - sign;
1620			return (multiply64((exp>-1) ? (((fh<<(1+exp))&0xFFFFFFFF)|((fl&0xFFFFFFFF)>>(31-exp))) : fh, 0xC90FDAA2)^sign) - sign;
1621		#endif
1622		}
1623
1624		/// Get arguments for atan2 function.
1625		/// \param abs half-precision floating-point value
1626		/// \return \a abs and sqrt(1 - \a abs^2) as Q0.30
1627		inline std::pair<uint32,uint32> atan2_args(unsigned int abs)
1628		{
1629			int exp = -15;
1630			for(; abs<0x400; abs<<=1,--exp) ;
1631			exp += abs >> 10;
1632			uint32 my = ((abs&0x3FF)|0x400) << 5, r = my * my;
1633			int rexp = 2 * exp;
1634			r = 0x40000000 - ((rexp>-31) ? ((r>>-rexp)|((r&((static_cast<uint32>(1)<<-rexp)-1))!=0)) : 1);
1635			for(rexp=0; r<0x40000000; r<<=1,--rexp) ;
1636			uint32 mx = sqrt<30>(r, rexp);
1637			int d = exp - rexp;
1638			if(d < 0)
1639				return std::make_pair((d<-14) ? ((my>>(-d-14))+((my>>(-d-15))&1)) : (my<<(14+d)), (mx<<14)+(r<<13)/mx);
1640			if(d > 0)
1641				return std::make_pair(my<<14, (d>14) ? ((mx>>(d-14))+((mx>>(d-15))&1)) : ((d==14) ? mx : ((mx<<(14-d))+(r<<(13-d))/mx)));
1642			return std::make_pair(my<<13, (mx<<13)+(r<<12)/mx);
1643		}
1644
1645		/// Get exponentials for hyperbolic computation
1646		/// \param abs half-precision floating-point value
1647		/// \param exp variable to take unbiased exponent of larger result
1648		/// \param n number of BKM iterations (at most 32)
1649		/// \return exp(abs) and exp(-\a abs) as Q1.31 with same exponent
1650		inline std::pair<uint32,uint32> hyperbolic_args(unsigned int abs, int &exp, unsigned int n = 32)
1651		{
1652			uint32 mx = detail::multiply64(static_cast<uint32>((abs&0x3FF)+((abs>0x3FF)<<10))<<21, 0xB8AA3B29), my;
1653			int e = (abs>>10) + (abs<=0x3FF);
1654			if(e < 14)
1655			{
1656				exp = 0;
1657				mx >>= 14 - e;
1658			}
1659			else
1660			{
1661				exp = mx >> (45-e);
1662				mx = (mx<<(e-14)) & 0x7FFFFFFF;
1663			}
1664			mx = exp2(mx, n);
1665			int d = exp << 1, s;
1666			if(mx > 0x80000000)
1667			{
1668				my = divide64(0x80000000, mx, s);
1669				my |= s;
1670				++d;
1671			}
1672			else
1673				my = mx;
1674			return std::make_pair(mx, (d<31) ? ((my>>d)|((my&((static_cast<uint32>(1)<<d)-1))!=0)) : 1);
1675		}
1676
1677		/// Postprocessing for binary exponential.
1678		/// \tparam R rounding mode to use
1679		/// \param m fractional part of as Q0.31
1680		/// \param exp absolute value of unbiased exponent
1681		/// \param esign sign of actual exponent
1682		/// \param sign sign bit of result
1683		/// \param n number of BKM iterations (at most 32)
1684		/// \return value converted to half-precision
1685		/// \exception FE_OVERFLOW on overflows
1686		/// \exception FE_UNDERFLOW on underflows
1687		/// \exception FE_INEXACT if value had to be rounded or \a I is `true`
1688		template<std::float_round_style R> unsigned int exp2_post(uint32 m, int exp, bool esign, unsigned int sign = 0, unsigned int n = 32)
1689		{
1690			if(esign)
1691			{
1692				exp = -exp - (m!=0);
1693				if(exp < -25)
1694					return underflow<R>(sign);
1695				else if(exp == -25)
1696					return rounded<R,false>(sign, 1, m!=0);
1697			}
1698			else if(exp > 15)
1699				return overflow<R>(sign);
1700			if(!m)
1701				return sign | (((exp+=15)>0) ? (exp<<10) : check_underflow(0x200>>-exp));
1702			m = exp2(m, n);
1703			int s = 0;
1704			if(esign)
1705				m = divide64(0x80000000, m, s);
1706			return fixed2half<R,31,false,false,true>(m, exp+14, sign, s);
1707		}
1708
1709		/// Postprocessing for binary logarithm.
1710		/// \tparam R rounding mode to use
1711		/// \tparam L logarithm for base transformation as Q1.31
1712		/// \param m fractional part of logarithm as Q0.31
1713		/// \param ilog signed integer part of logarithm
1714		/// \param exp biased exponent of result
1715		/// \param sign sign bit of result
1716		/// \return value base-transformed and converted to half-precision
1717		/// \exception FE_OVERFLOW on overflows
1718		/// \exception FE_UNDERFLOW on underflows
1719		/// \exception FE_INEXACT if no other exception occurred
1720		template<std::float_round_style R,uint32 L> unsigned int log2_post(uint32 m, int ilog, int exp, unsigned int sign = 0)
1721		{
1722			uint32 msign = sign_mask(ilog);
1723			m = (((static_cast<uint32>(ilog)<<27)+(m>>4))^msign) - msign;
1724			if(!m)
1725				return 0;
1726			for(; m<0x80000000; m<<=1,--exp) ;
1727			int i = m >= L, s;
1728			exp += i;
1729			m >>= 1 + i;
1730			sign ^= msign & 0x8000;
1731			if(exp < -11)
1732				return underflow<R>(sign);
1733			m = divide64(m, L, s);
1734			return fixed2half<R,30,false,false,true>(m, exp, sign, 1);
1735		}
1736
1737		/// Hypotenuse square root and postprocessing.
1738		/// \tparam R rounding mode to use
1739		/// \param r mantissa as Q2.30
1740		/// \param exp biased exponent
1741		/// \return square root converted to half-precision
1742		/// \exception FE_OVERFLOW on overflows
1743		/// \exception FE_UNDERFLOW on underflows
1744		/// \exception FE_INEXACT if value had to be rounded
1745		template<std::float_round_style R> unsigned int hypot_post(uint32 r, int exp)
1746		{
1747			int i = r >> 31;
1748			if((exp+=i) > 46)
1749				return overflow<R>();
1750			if(exp < -34)
1751				return underflow<R>();
1752			r = (r>>i) | (r&i);
1753			uint32 m = sqrt<30>(r, exp+=15);
1754			return fixed2half<R,15,false,false,false>(m, exp-1, 0, r!=0);
1755		}
1756
1757		/// Division and postprocessing for tangents.
1758		/// \tparam R rounding mode to use
1759		/// \param my dividend as Q1.31
1760		/// \param mx divisor as Q1.31
1761		/// \param exp biased exponent of result
1762		/// \param sign sign bit of result
1763		/// \return quotient converted to half-precision
1764		/// \exception FE_OVERFLOW on overflows
1765		/// \exception FE_UNDERFLOW on underflows
1766		/// \exception FE_INEXACT if no other exception occurred
1767		template<std::float_round_style R> unsigned int tangent_post(uint32 my, uint32 mx, int exp, unsigned int sign = 0)
1768		{
1769			int i = my >= mx, s;
1770			exp += i;
1771			if(exp > 29)
1772				return overflow<R>(sign);
1773			if(exp < -11)
1774				return underflow<R>(sign);
1775			uint32 m = divide64(my>>(i+1), mx, s);
1776			return fixed2half<R,30,false,false,true>(m, exp, sign, s);
1777		}
1778
1779		/// Area function and postprocessing.
1780		/// This computes the value directly in Q2.30 using the representation `asinh|acosh(x) = log(x+sqrt(x^2+|-1))`.
1781		/// \tparam R rounding mode to use
1782		/// \tparam S `true` for asinh, `false` for acosh
1783		/// \param arg half-precision argument
1784		/// \return asinh|acosh(\a arg) converted to half-precision
1785		/// \exception FE_OVERFLOW on overflows
1786		/// \exception FE_UNDERFLOW on underflows
1787		/// \exception FE_INEXACT if no other exception occurred
1788		template<std::float_round_style R,bool S> unsigned int area(unsigned int arg)
1789		{
1790			int abs = arg & 0x7FFF, expx = (abs>>10) + (abs<=0x3FF) - 15, expy = -15, ilog, i;
1791			uint32 mx = static_cast<uint32>((abs&0x3FF)|((abs>0x3FF)<<10)) << 20, my, r;
1792			for(; abs<0x400; abs<<=1,--expy) ;
1793			expy += abs >> 10;
1794			r = ((abs&0x3FF)|0x400) << 5;
1795			r *= r;
1796			i = r >> 31;
1797			expy = 2*expy + i;
1798			r >>= i;
1799			if(S)
1800			{
1801				if(expy < 0)
1802				{
1803					r = 0x40000000 + ((expy>-30) ? ((r>>-expy)|((r&((static_cast<uint32>(1)<<-expy)-1))!=0)) : 1);
1804					expy = 0;
1805				}
1806				else
1807				{
1808					r += 0x40000000 >> expy;
1809					i = r >> 31;
1810					r = (r>>i) | (r&i);
1811					expy += i;
1812				}
1813			}
1814			else
1815			{
1816				r -= 0x40000000 >> expy;
1817				for(; r<0x40000000; r<<=1,--expy) ;
1818			}
1819			my = sqrt<30>(r, expy);
1820			my = (my<<15) + (r<<14)/my;
1821			if(S)
1822			{
1823				mx >>= expy - expx;
1824				ilog = expy;
1825			}
1826			else
1827			{
1828				my >>= expx - expy;
1829				ilog = expx;
1830			}
1831			my += mx;
1832			i = my >> 31;
1833			static const int G = S && (R==std::round_to_nearest);
1834			return log2_post<R,0xB8AA3B2A>(log2(my>>i, 26+S+G)+(G<<3), ilog+i, 17, arg&(static_cast<unsigned>(S)<<15));
1835		}
1836
1837		/// Class for 1.31 unsigned floating-point computation
1838		struct f31
1839		{
1840			/// Constructor.
1841			/// \param mant mantissa as 1.31
1842			/// \param e exponent
1843			HALF_CONSTEXPR f31(uint32 mant, int e) : m(mant), exp(e) {}
1844
1845			/// Constructor.
1846			/// \param abs unsigned half-precision value
1847			f31(unsigned int abs) : exp(-15)
1848			{
1849				for(; abs<0x400; abs<<=1,--exp) ;
1850				m = static_cast<uint32>((abs&0x3FF)|0x400) << 21;
1851				exp += (abs>>10);
1852			}
1853
1854			/// Addition operator.
1855			/// \param a first operand
1856			/// \param b second operand
1857			/// \return \a a + \a b
1858			friend f31 operator+(f31 a, f31 b)
1859			{
1860				if(b.exp > a.exp)
1861					std::swap(a, b);
1862				int d = a.exp - b.exp;
1863				uint32 m = a.m + ((d<32) ? (b.m>>d) : 0);
1864				int i = (m&0xFFFFFFFF) < a.m;
1865				return f31(((m+i)>>i)|0x80000000, a.exp+i);
1866			}
1867
1868			/// Subtraction operator.
1869			/// \param a first operand
1870			/// \param b second operand
1871			/// \return \a a - \a b
1872			friend f31 operator-(f31 a, f31 b)
1873			{
1874				int d = a.exp - b.exp, exp = a.exp;
1875				uint32 m = a.m - ((d<32) ? (b.m>>d) : 0);
1876				if(!m)
1877					return f31(0, -32);
1878				for(; m<0x80000000; m<<=1,--exp) ;
1879				return f31(m, exp);
1880			}
1881
1882			/// Multiplication operator.
1883			/// \param a first operand
1884			/// \param b second operand
1885			/// \return \a a * \a b
1886			friend f31 operator*(f31 a, f31 b)
1887			{
1888				uint32 m = multiply64(a.m, b.m);
1889				int i = m >> 31;
1890				return f31(m<<(1-i), a.exp + b.exp + i);
1891			}
1892
1893			/// Division operator.
1894			/// \param a first operand
1895			/// \param b second operand
1896			/// \return \a a / \a b
1897			friend f31 operator/(f31 a, f31 b)
1898			{
1899				int i = a.m >= b.m, s;
1900				uint32 m = divide64((a.m+i)>>i, b.m, s);
1901				return f31(m, a.exp - b.exp + i - 1);
1902			}
1903
1904			uint32 m;			///< mantissa as 1.31.
1905			int exp;			///< exponent.
1906		};
1907
1908		/// Error function and postprocessing.
1909		/// This computes the value directly in Q1.31 using the approximations given 
1910		/// [here](https://en.wikipedia.org/wiki/Error_function#Approximation_with_elementary_functions).
1911		/// \tparam R rounding mode to use
1912		/// \tparam C `true` for comlementary error function, `false` else
1913		/// \param arg half-precision function argument
1914		/// \return approximated value of error function in half-precision
1915		/// \exception FE_OVERFLOW on overflows
1916		/// \exception FE_UNDERFLOW on underflows
1917		/// \exception FE_INEXACT if no other exception occurred
1918		template<std::float_round_style R,bool C> unsigned int erf(unsigned int arg)
1919		{
1920			unsigned int abs = arg & 0x7FFF, sign = arg & 0x8000;
1921			f31 x(abs), x2 = x * x * f31(0xB8AA3B29, 0), t = f31(0x80000000, 0) / (f31(0x80000000, 0)+f31(0xA7BA054A, -2)*x), t2 = t * t;
1922			f31 e = ((f31(0x87DC2213, 0)*t2+f31(0xB5F0E2AE, 0))*t2+f31(0x82790637, -2)-(f31(0xBA00E2B8, 0)*t2+f31(0x91A98E62, -2))*t) * t /
1923					((x2.exp<0) ? f31(exp2((x2.exp>-32) ? (x2.m>>-x2.exp) : 0, 30), 0) : f31(exp2((x2.m<<x2.exp)&0x7FFFFFFF, 22), x2.m>>(31-x2.exp)));
1924			return (!C || sign) ? fixed2half<R,31,false,true,true>(0x80000000-(e.m>>(C-e.exp)), 14+C, sign&(C-1U)) :
1925					(e.exp<-25) ? underflow<R>() : fixed2half<R,30,false,false,true>(e.m>>1, e.exp+14, 0, e.m&1);
1926		}
1927
1928		/// Gamma function and postprocessing.
1929		/// This approximates the value of either the gamma function or its logarithm directly in Q1.31.
1930		/// \tparam R rounding mode to use
1931		/// \tparam L `true` for lograithm of gamma function, `false` for gamma function
1932		/// \param arg half-precision floating-point value
1933		/// \return lgamma/tgamma(\a arg) in half-precision
1934		/// \exception FE_OVERFLOW on overflows
1935		/// \exception FE_UNDERFLOW on underflows
1936		/// \exception FE_INEXACT if \a arg is not a positive integer
1937		template<std::float_round_style R,bool L> unsigned int gamma(unsigned int arg)
1938		{
1939/*			static const double p[] ={ 2.50662827563479526904, 225.525584619175212544, -268.295973841304927459, 80.9030806934622512966, -5.00757863970517583837, 0.0114684895434781459556 };
1940			double t = arg + 4.65, s = p[0];
1941			for(unsigned int i=0; i<5; ++i)
1942				s += p[i+1] / (arg+i);
1943			return std::log(s) + (arg-0.5)*std::log(t) - t;
1944*/			static const f31 pi(0xC90FDAA2, 1), lbe(0xB8AA3B29, 0);
1945			unsigned int abs = arg & 0x7FFF, sign = arg & 0x8000;
1946			bool bsign = sign != 0;
1947			f31 z(abs), x = sign ? (z+f31(0x80000000, 0)) : z, t = x + f31(0x94CCCCCD, 2), s =
1948				f31(0xA06C9901, 1) + f31(0xBBE654E2, -7)/(x+f31(0x80000000, 2)) + f31(0xA1CE6098, 6)/(x+f31(0x80000000, 1))
1949				+ f31(0xE1868CB7, 7)/x - f31(0x8625E279, 8)/(x+f31(0x80000000, 0)) - f31(0xA03E158F, 2)/(x+f31(0xC0000000, 1));
1950			int i = (s.exp>=2) + (s.exp>=4) + (s.exp>=8) + (s.exp>=16);
1951			s = f31((static_cast<uint32>(s.exp)<<(31-i))+(log2(s.m>>1, 28)>>i), i) / lbe;
1952			if(x.exp != -1 || x.m != 0x80000000)
1953			{
1954				i = (t.exp>=2) + (t.exp>=4) + (t.exp>=8);
1955				f31 l = f31((static_cast<uint32>(t.exp)<<(31-i))+(log2(t.m>>1, 30)>>i), i) / lbe;
1956				s = (x.exp<-1) ? (s-(f31(0x80000000, -1)-x)*l) : (s+(x-f31(0x80000000, -1))*l);
1957			}
1958			s = x.exp ? (s-t) : (t-s);
1959			if(bsign)
1960			{
1961				if(z.exp >= 0)
1962				{
1963					sign &= (L|((z.m>>(31-z.exp))&1)) - 1;
1964					for(z=f31((z.m<<(1+z.exp))&0xFFFFFFFF, -1); z.m<0x80000000; z.m<<=1,--z.exp) ;
1965				}
1966				if(z.exp == -1)
1967					z = f31(0x80000000, 0) - z;
1968				if(z.exp < -1)
1969				{
1970					z = z * pi;
1971					z.m = sincos(z.m>>(1-z.exp), 30).first;
1972					for(z.exp=1; z.m<0x80000000; z.m<<=1,--z.exp) ;
1973				}
1974				else
1975					z = f31(0x80000000, 0);
1976			}
1977			if(L)
1978			{
1979				if(bsign)
1980				{
1981					f31 l(0x92868247, 0);
1982					if(z.exp < 0)
1983					{
1984						uint32 m = log2((z.m+1)>>1, 27);
1985						z = f31(-((static_cast<uint32>(z.exp)<<26)+(m>>5)), 5);
1986						for(; z.m<0x80000000; z.m<<=1,--z.exp) ;
1987						l = l + z / lbe;
1988					}
1989					sign = static_cast<unsigned>(x.exp&&(l.exp<s.exp||(l.exp==s.exp&&l.m<s.m))) << 15;
1990					s = sign ? (s-l) : x.exp ? (l-s) : (l+s);
1991				}
1992				else
1993				{
1994					sign = static_cast<unsigned>(x.exp==0) << 15;
1995					if(s.exp < -24)
1996						return underflow<R>(sign);
1997					if(s.exp > 15)
1998						return overflow<R>(sign);
1999				}
2000			}
2001			else
2002			{
2003				s = s * lbe;
2004				uint32 m;
2005				if(s.exp < 0)
2006				{
2007					m = s.m >> -s.exp;
2008					s.exp = 0;
2009				}
2010				else
2011				{
2012					m = (s.m<<s.exp) & 0x7FFFFFFF;
2013					s.exp = (s.m>>(31-s.exp));
2014				}
2015				s.m = exp2(m, 27);
2016				if(!x.exp)
2017					s = f31(0x80000000, 0) / s;
2018				if(bsign)
2019				{
2020					if(z.exp < 0)
2021						s = s * z;
2022					s = pi / s;
2023					if(s.exp < -24)
2024						return underflow<R>(sign);
2025				}
2026				else if(z.exp > 0 && !(z.m&((1<<(31-z.exp))-1)))
2027					return ((s.exp+14)<<10) + (s.m>>21);
2028				if(s.exp > 15)
2029					return overflow<R>(sign);
2030			}
2031			return fixed2half<R,31,false,false,true>(s.m, s.exp+14, sign);
2032		}
2033		/// \}
2034
2035		template<typename,typename,std::float_round_style> struct half_caster;
2036	}
2037
2038	/// Half-precision floating-point type.
2039	/// This class implements an IEEE-conformant half-precision floating-point type with the usual arithmetic 
2040	/// operators and conversions. It is implicitly convertible to single-precision floating-point, which makes artihmetic 
2041	/// expressions and functions with mixed-type operands to be of the most precise operand type.
2042	///
2043	/// According to the C++98/03 definition, the half type is not a POD type. But according to C++11's less strict and 
2044	/// extended definitions it is both a standard layout type and a trivially copyable type (even if not a POD type), which 
2045	/// means it can be standard-conformantly copied using raw binary copies. But in this context some more words about the 
2046	/// actual size of the type. Although the half is representing an IEEE 16-bit type, it does not neccessarily have to be of 
2047	/// exactly 16-bits size. But on any reasonable implementation the actual binary representation of this type will most 
2048	/// probably not ivolve any additional "magic" or padding beyond the simple binary representation of the underlying 16-bit 
2049	/// IEEE number, even if not strictly guaranteed by the standard. But even then it only has an actual size of 16 bits if 
2050	/// your C++ implementation supports an unsigned integer type of exactly 16 bits width. But this should be the case on 
2051	/// nearly any reasonable platform.
2052	///
2053	/// So if your C++ implementation is not totally exotic or imposes special alignment requirements, it is a reasonable 
2054	/// assumption that the data of a half is just comprised of the 2 bytes of the underlying IEEE representation.
2055	class half
2056	{
2057	public:
2058		/// \name Construction and assignment
2059		/// \{
2060
2061		/// Default constructor.
2062		/// This initializes the half to 0. Although this does not match the builtin types' default-initialization semantics 
2063		/// and may be less efficient than no initialization, it is needed to provide proper value-initialization semantics.
2064		HALF_CONSTEXPR half() HALF_NOEXCEPT : data_() {}
2065
2066		/// Conversion constructor.
2067		/// \param rhs float to convert
2068		/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
2069		explicit half(float rhs) : data_(static_cast<detail::uint16>(detail::float2half<round_style>(rhs))) {}
2070	
2071		/// Conversion to single-precision.
2072		/// \return single precision value representing expression value
2073		operator float() const { return detail::half2float<float>(data_); }
2074
2075		/// Assignment operator.
2076		/// \param rhs single-precision value to copy from
2077		/// \return reference to this half
2078		/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
2079		half& operator=(float rhs) { data_ = static_cast<detail::uint16>(detail::float2half<round_style>(rhs)); return *this; }
2080
2081		/// \}
2082		/// \name Arithmetic updates
2083		/// \{
2084
2085		/// Arithmetic assignment.
2086		/// \tparam T type of concrete half expression
2087		/// \param rhs half expression to add
2088		/// \return reference to this half
2089		/// \exception FE_... according to operator+(half,half)
2090		half& operator+=(half rhs) { return *this = *this + rhs; }
2091
2092		/// Arithmetic assignment.
2093		/// \tparam T type of concrete half expression
2094		/// \param rhs half expression to subtract
2095		/// \return reference to this half
2096		/// \exception FE_... according to operator-(half,half)
2097		half& operator-=(half rhs) { return *this = *this - rhs; }
2098
2099		/// Arithmetic assignment.
2100		/// \tparam T type of concrete half expression
2101		/// \param rhs half expression to multiply with
2102		/// \return reference to this half
2103		/// \exception FE_... according to operator*(half,half)
2104		half& operator*=(half rhs) { return *this = *this * rhs; }
2105
2106		/// Arithmetic assignment.
2107		/// \tparam T type of concrete half expression
2108		/// \param rhs half expression to divide by
2109		/// \return reference to this half
2110		/// \exception FE_... according to operator/(half,half)
2111		half& operator/=(half rhs) { return *this = *this / rhs; }
2112
2113		/// Arithmetic assignment.
2114		/// \param rhs single-precision value to add
2115		/// \return reference to this half
2116		/// \exception FE_... according to operator=()
2117		half& operator+=(float rhs) { return *this = *this + rhs; }
2118
2119		/// Arithmetic assignment.
2120		/// \param rhs single-precision value to subtract
2121		/// \return reference to this half
2122		/// \exception FE_... according to operator=()
2123		half& operator-=(float rhs) { return *this = *this - rhs; }
2124
2125		/// Arithmetic assignment.
2126		/// \param rhs single-precision value to multiply with
2127		/// \return reference to this half
2128		/// \exception FE_... according to operator=()
2129		half& operator*=(float rhs) { return *this = *this * rhs; }
2130
2131		/// Arithmetic assignment.
2132		/// \param rhs single-precision value to divide by
2133		/// \return reference to this half
2134		/// \exception FE_... according to operator=()
2135		half& operator/=(float rhs) { return *this = *this / rhs; }
2136
2137		/// \}
2138		/// \name Increment and decrement
2139		/// \{
2140
2141		/// Prefix increment.
2142		/// \return incremented half value
2143		/// \exception FE_... according to operator+(half,half)
2144		half& operator++() { return *this = *this + half(detail::binary, 0x3C00); }
2145
2146		/// Prefix decrement.
2147		/// \return decremented half value
2148		/// \exception FE_... according to operator-(half,half)
2149		half& operator--() { return *this = *this + half(detail::binary, 0xBC00); }
2150
2151		/// Postfix increment.
2152		/// \return non-incremented half value
2153		/// \exception FE_... according to operator+(half,half)
2154		half operator++(int) { half out(*this); ++*this; return out; }
2155
2156		/// Postfix decrement.
2157		/// \return non-decremented half value
2158		/// \exception FE_... according to operator-(half,half)
2159		half operator--(int) { half out(*this); --*this; return out; }
2160		/// \}
2161	
2162	private:
2163		/// Rounding mode to use
2164		static const std::float_round_style round_style = (std::float_round_style)(HALF_ROUND_STYLE);
2165
2166		/// Constructor.
2167		/// \param bits binary representation to set half to
2168		HALF_CONSTEXPR half(detail::binary_t, unsigned int bits) HALF_NOEXCEPT : data_(static_cast<detail::uint16>(bits)) {}
2169
2170		/// Internal binary representation
2171		detail::uint16 data_;
2172
2173	#ifndef HALF_DOXYGEN_ONLY
2174		friend HALF_CONSTEXPR_NOERR bool operator==(half, half);
2175		friend HALF_CONSTEXPR_NOERR bool operator!=(half, half);
2176		friend HALF_CONSTEXPR_NOERR bool operator<(half, half);
2177		friend HALF_CONSTEXPR_NOERR bool operator>(half, half);
2178		friend HALF_CONSTEXPR_NOERR bool operator<=(half, half);
2179		friend HALF_CONSTEXPR_NOERR bool operator>=(half, half);
2180		friend HALF_CONSTEXPR half operator-(half);
2181		friend half operator+(half, half);
2182		friend half operator-(half, half);
2183		friend half operator*(half, half);
2184		friend half operator/(half, half);
2185		template<typename charT,typename traits> friend std::basic_ostream<charT,traits>& operator<<(std::basic_ostream<charT,traits>&, half);
2186		template<typename charT,typename traits> friend std::basic_istream<charT,traits>& operator>>(std::basic_istream<charT,traits>&, half&);
2187		friend HALF_CONSTEXPR half fabs(half);
2188		friend half fmod(half, half);
2189		friend half remainder(half, half);
2190		friend half remquo(half, half, int*);
2191		friend half fma(half, half, half);
2192		friend HALF_CONSTEXPR_NOERR half fmax(half, half);
2193		friend HALF_CONSTEXPR_NOERR half fmin(half, half);
2194		friend half fdim(half, half);
2195		friend half nanh(const char*);
2196		friend half exp(half);
2197		friend half exp2(half);
2198		friend half expm1(half);
2199		friend half log(half);
2200		friend half log10(half);
2201		friend half log2(half);
2202		friend half log1p(half);
2203		friend half sqrt(half);
2204		friend half rsqrt(half);
2205		friend half cbrt(half);
2206		friend half hypot(half, half);
2207		friend half hypot(half, half, half);
2208		friend half pow(half, half);
2209		friend void sincos(half, half*, half*);
2210		friend half sin(half);
2211		friend half cos(half);
2212		friend half tan(half);
2213		friend half asin(half);
2214		friend half acos(half);
2215		friend half atan(half);
2216		friend half atan2(half, half);
2217		friend half sinh(half);
2218		friend half cosh(half);
2219		friend half tanh(half);
2220		friend half asinh(half);
2221		friend half acosh(half);
2222		friend half atanh(half);
2223		friend half erf(half);
2224		friend half erfc(half);
2225		friend half lgamma(half);
2226		friend half tgamma(half);
2227		friend half ceil(half);
2228		friend half floor(half);
2229		friend half trunc(half);
2230		friend half round(half);
2231		friend long lround(half);
2232		friend half rint(half);
2233		friend long lrint(half);
2234		friend half nearbyint(half);
2235	#ifdef HALF_ENABLE_CPP11_LONG_LONG
2236		friend long long llround(half);
2237		friend long long llrint(half);
2238	#endif
2239		friend half frexp(half, int*);
2240		friend half scalbln(half, long);
2241		friend half modf(half, half*);
2242		friend int ilogb(half);
2243		friend half logb(half);
2244		friend half nextafter(half, half);
2245		friend half nexttoward(half, long double);
2246		friend HALF_CONSTEXPR half copysign(half, half);
2247		friend HALF_CONSTEXPR int fpclassify(half);
2248		friend HALF_CONSTEXPR bool isfinite(half);
2249		friend HALF_CONSTEXPR bool isinf(half);
2250		friend HALF_CONSTEXPR bool isnan(half);
2251		friend HALF_CONSTEXPR bool isnormal(half);
2252		friend HALF_CONSTEXPR bool signbit(half);
2253		friend HALF_CONSTEXPR bool isgreater(half, half);
2254		friend HALF_CONSTEXPR bool isgreaterequal(half, half);
2255		friend HALF_CONSTEXPR bool isless(half, half);
2256		friend HALF_CONSTEXPR bool islessequal(half, half);
2257		friend HALF_CONSTEXPR bool islessgreater(half, half);
2258		template<typename,typename,std::float_round_style> friend struct detail::half_caster;
2259		friend class std::numeric_limits<half>;
2260	#if HALF_ENABLE_CPP11_HASH
2261		friend struct std::hash<half>;
2262	#endif
2263	#if HALF_ENABLE_CPP11_USER_LITERALS
2264		friend half literal::operator "" _h(long double);
2265	#endif
2266	#endif
2267	};
2268
2269#if HALF_ENABLE_CPP11_USER_LITERALS
2270	namespace literal
2271	{
2272		/// Half literal.
2273		/// While this returns a properly rounded half-precision value, half literals can unfortunately not be constant 
2274		/// expressions due to rather involved conversions. So don't expect this to be a literal literal without involving 
2275		/// conversion operations at runtime. It is a convenience feature, not a performance optimization.
2276		/// \param value literal value
2277		/// \return half with of given value (possibly rounded)
2278		/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
2279		inline half operator "" _h(long double value) { return half(detail::binary, detail::float2half<half::round_style>(value)); }
2280	}
2281#endif
2282
2283	namespace detail
2284	{
2285		/// Helper class for half casts.
2286		/// This class template has to be specialized for all valid cast arguments to define an appropriate static 
2287		/// `cast` member function and a corresponding `type` member denoting its return type.
2288		/// \tparam T destination type
2289		/// \tparam U source type
2290		/// \tparam R rounding mode to use
2291		template<typename T,typename U,std::float_round_style R=(std::float_round_style)(HALF_ROUND_STYLE)> struct half_caster {};
2292		template<typename U,std::float_round_style R> struct half_caster<half,U,R>
2293		{
2294		#if HALF_ENABLE_CPP11_STATIC_ASSERT && HALF_ENABLE_CPP11_TYPE_TRAITS
2295			static_assert(std::is_arithmetic<U>::value, "half_cast from non-arithmetic type unsupported");
2296		#endif
2297
2298			static half cast(U arg) { return cast_impl(arg, is_float<U>()); };
2299
2300		private:
2301			static half cast_impl(U arg, true_type) { return half(binary, float2half<R>(arg)); }
2302			static half cast_impl(U arg, false_type) { return half(binary, int2half<R>(arg)); }
2303		};
2304		template<typename T,std::float_round_style R> struct half_caster<T,half,R>
2305		{
2306		#if HALF_ENABLE_CPP11_STATIC_ASSERT && HALF_ENABLE_CPP11_TYPE_TRAITS
2307			static_assert(std::is_arithmetic<T>::value, "half_cast to non-arithmetic type unsupported");
2308		#endif
2309
2310			static T cast(half arg) { return cast_impl(arg, is_float<T>()); }
2311
2312		private:
2313			static T cast_impl(half arg, true_type) { return half2float<T>(arg.data_); }
2314			static T cast_impl(half arg, false_type) { return half2int<R,true,true,T>(arg.data_); }
2315		};
2316		template<std::float_round_style R> struct half_caster<half,half,R>
2317		{
2318			static half cast(half arg) { return arg; }
2319		};
2320	}
2321}
2322
2323/// Extensions to the C++ standard library.
2324namespace std
2325{
2326	/// Numeric limits for half-precision floats.
2327	/// **See also:** Documentation for [std::numeric_limits](https://en.cppreference.com/w/cpp/types/numeric_limits)
2328	template<> class numeric_limits<half_float::half>
2329	{
2330	public:
2331		/// Is template specialization.
2332		static HALF_CONSTEXPR_CONST bool is_specialized = true;
2333
2334		/// Supports signed values.
2335		static HALF_CONSTEXPR_CONST bool is_signed = true;
2336
2337		/// Is not an integer type.
2338		static HALF_CONSTEXPR_CONST bool is_integer = false;
2339
2340		/// Is not exact.
2341		static HALF_CONSTEXPR_CONST bool is_exact = false;
2342
2343		/// Doesn't provide modulo arithmetic.
2344		static HALF_CONSTEXPR_CONST bool is_modulo = false;
2345
2346		/// Has a finite set of values.
2347		static HALF_CONSTEXPR_CONST bool is_bounded = true;
2348
2349		/// IEEE conformant.
2350		static HALF_CONSTEXPR_CONST bool is_iec559 = true;
2351
2352		/// Supports infinity.
2353		static HALF_CONSTEXPR_CONST bool has_infinity = true;
2354
2355		/// Supports quiet NaNs.
2356		static HALF_CONSTEXPR_CONST bool has_quiet_NaN = true;
2357
2358		/// Supports signaling NaNs.
2359		static HALF_CONSTEXPR_CONST bool has_signaling_NaN = true;
2360
2361		/// Supports subnormal values.
2362		static HALF_CONSTEXPR_CONST float_denorm_style has_denorm = denorm_present;
2363
2364		/// Supports no denormalization detection.
2365		static HALF_CONSTEXPR_CONST bool has_denorm_loss = false;
2366
2367	#if HALF_ERRHANDLING_THROWS
2368		static HALF_CONSTEXPR_CONST bool traps = true;
2369	#else
2370		/// Traps only if [HALF_ERRHANDLING_THROW_...](\ref HALF_ERRHANDLING_THROW_INVALID) is acitvated.
2371		static HALF_CONSTEXPR_CONST bool traps = false;
2372	#endif
2373
2374		/// Does not support no pre-rounding underflow detection.
2375		static HALF_CONSTEXPR_CONST bool tinyness_before = false;
2376
2377		/// Rounding mode.
2378		static HALF_CONSTEXPR_CONST float_round_style round_style = half_float::half::round_style;
2379
2380		/// Significant digits.
2381		static HALF_CONSTEXPR_CONST int digits = 11;
2382
2383		/// Significant decimal digits.
2384		static HALF_CONSTEXPR_CONST int digits10 = 3;
2385
2386		/// Required decimal digits to represent all possible values.
2387		static HALF_CONSTEXPR_CONST int max_digits10 = 5;
2388
2389		/// Number base.
2390		static HALF_CONSTEXPR_CONST int radix = 2;
2391
2392		/// One more than smallest exponent.
2393		static HALF_CONSTEXPR_CONST int min_exponent = -13;
2394
2395		/// Smallest normalized representable power of 10.
2396		static HALF_CONSTEXPR_CONST int min_exponent10 = -4;
2397
2398		/// One more than largest exponent
2399		static HALF_CONSTEXPR_CONST int max_exponent = 16;
2400
2401		/// Largest finitely representable power of 10.
2402		static HALF_CONSTEXPR_CONST int max_exponent10 = 4;
2403
2404		/// Smallest positive normal value.
2405		static HALF_CONSTEXPR half_float::half min() HALF_NOTHROW { return half_float::half(half_float::detail::binary, 0x0400); }
2406
2407		/// Smallest finite value.
2408		static HALF_CONSTEXPR half_float::half lowest() HALF_NOTHROW { return half_float::half(half_float::detail::binary, 0xFBFF); }
2409
2410		/// Largest finite value.
2411		static HALF_CONSTEXPR half_float::half max() HALF_NOTHROW { return half_float::half(half_float::detail::binary, 0x7BFF); }
2412
2413		/// Difference between 1 and next representable value.
2414		static HALF_CONSTEXPR half_float::half epsilon() HALF_NOTHROW { return half_float::half(half_float::detail::binary, 0x1400); }
2415
2416		/// Maximum rounding error in ULP (units in the last place).
2417		static HALF_CONSTEXPR half_float::half round_error() HALF_NOTHROW
2418			{ return half_float::half(half_float::detail::binary, (round_style==std::round_to_nearest) ? 0x3800 : 0x3C00); }
2419
2420		/// Positive infinity.
2421		static HALF_CONSTEXPR half_float::half infinity() HALF_NOTHROW { return half_float::half(half_float::detail::binary, 0x7C00); }
2422
2423		/// Quiet NaN.
2424		static HALF_CONSTEXPR half_float::half quiet_NaN() HALF_NOTHROW { return half_float::half(half_float::detail::binary, 0x7FFF); }
2425
2426		/// Signaling NaN.
2427		static HALF_CONSTEXPR half_float::half signaling_NaN() HALF_NOTHROW { return half_float::half(half_float::detail::binary, 0x7DFF); }
2428
2429		/// Smallest positive subnormal value.
2430		static HALF_CONSTEXPR half_float::half denorm_min() HALF_NOTHROW { return half_float::half(half_float::detail::binary, 0x0001); }
2431	};
2432
2433#if HALF_ENABLE_CPP11_HASH
2434	/// Hash function for half-precision floats.
2435	/// This is only defined if C++11 `std::hash` is supported and enabled.
2436	///
2437	/// **See also:** Documentation for [std::hash](https://en.cppreference.com/w/cpp/utility/hash)
2438	template<> struct hash<half_float::half>
2439	{
2440		/// Type of function argument.
2441		typedef half_float::half argument_type;
2442
2443		/// Function return type.
2444		typedef size_t result_type;
2445
2446		/// Compute hash function.
2447		/// \param arg half to hash
2448		/// \return hash value
2449		result_type operator()(argument_type arg) const { return hash<half_float::detail::uint16>()(arg.data_&-static_cast<unsigned>(arg.data_!=0x8000)); }
2450	};
2451#endif
2452}
2453
2454namespace half_float
2455{
2456	/// \anchor compop
2457	/// \name Comparison operators
2458	/// \{
2459
2460	/// Comparison for equality.
2461	/// \param x first operand
2462	/// \param y second operand
2463	/// \retval true if operands equal
2464	/// \retval false else
2465	/// \exception FE_INVALID if \a x or \a y is NaN
2466	inline HALF_CONSTEXPR_NOERR bool operator==(half x, half y)
2467	{
2468		return !detail::compsignal(x.data_, y.data_) && (x.data_==y.data_ || !((x.data_|y.data_)&0x7FFF));
2469	}
2470
2471	/// Comparison for inequality.
2472	/// \param x first operand
2473	/// \param y second operand
2474	/// \retval true if operands not equal
2475	/// \retval false else
2476	/// \exception FE_INVALID if \a x or \a y is NaN
2477	inline HALF_CONSTEXPR_NOERR bool operator!=(half x, half y)
2478	{
2479		return detail::compsignal(x.data_, y.data_) || (x.data_!=y.data_ && ((x.data_|y.data_)&0x7FFF));
2480	}
2481
2482	/// Comparison for less than.
2483	/// \param x first operand
2484	/// \param y second operand
2485	/// \retval true if \a x less than \a y
2486	/// \retval false else
2487	/// \exception FE_INVALID if \a x or \a y is NaN
2488	inline HALF_CONSTEXPR_NOERR bool operator<(half x, half y)
2489	{
2490		return !detail::compsignal(x.data_, y.data_) &&
2491			((x.data_^(0x8000|(0x8000-(x.data_>>15))))+(x.data_>>15)) < ((y.data_^(0x8000|(0x8000-(y.data_>>15))))+(y.data_>>15));
2492	}
2493
2494	/// Comparison for greater than.
2495	/// \param x first operand
2496	/// \param y second operand
2497	/// \retval true if \a x greater than \a y
2498	/// \retval false else
2499	/// \exception FE_INVALID if \a x or \a y is NaN
2500	inline HALF_CONSTEXPR_NOERR bool operator>(half x, half y)
2501	{
2502		return !detail::compsignal(x.data_, y.data_) &&
2503			((x.data_^(0x8000|(0x8000-(x.data_>>15))))+(x.data_>>15)) > ((y.data_^(0x8000|(0x8000-(y.data_>>15))))+(y.data_>>15));
2504	}
2505
2506	/// Comparison for less equal.
2507	/// \param x first operand
2508	/// \param y second operand
2509	/// \retval true if \a x less equal \a y
2510	/// \retval false else
2511	/// \exception FE_INVALID if \a x or \a y is NaN
2512	inline HALF_CONSTEXPR_NOERR bool operator<=(half x, half y)
2513	{
2514		return !detail::compsignal(x.data_, y.data_) &&
2515			((x.data_^(0x8000|(0x8000-(x.data_>>15))))+(x.data_>>15)) <= ((y.data_^(0x8000|(0x8000-(y.data_>>15))))+(y.data_>>15));
2516	}
2517
2518	/// Comparison for greater equal.
2519	/// \param x first operand
2520	/// \param y second operand
2521	/// \retval true if \a x greater equal \a y
2522	/// \retval false else
2523	/// \exception FE_INVALID if \a x or \a y is NaN
2524	inline HALF_CONSTEXPR_NOERR bool operator>=(half x, half y)
2525	{
2526		return !detail::compsignal(x.data_, y.data_) &&
2527			((x.data_^(0x8000|(0x8000-(x.data_>>15))))+(x.data_>>15)) >= ((y.data_^(0x8000|(0x8000-(y.data_>>15))))+(y.data_>>15));
2528	}
2529
2530	/// \}
2531	/// \anchor arithmetics
2532	/// \name Arithmetic operators
2533	/// \{
2534
2535	/// Identity.
2536	/// \param arg operand
2537	/// \return unchanged operand
2538	inline HALF_CONSTEXPR half operator+(half arg) { return arg; }
2539
2540	/// Negation.
2541	/// \param arg operand
2542	/// \return negated operand
2543	inline HALF_CONSTEXPR half operator-(half arg) { return half(detail::binary, arg.data_^0x8000); }
2544
2545	/// Addition.
2546	/// This operation is exact to rounding for all rounding modes.
2547	/// \param x left operand
2548	/// \param y right operand
2549	/// \return sum of half expressions
2550	/// \exception FE_INVALID if \a x and \a y are infinities with different signs or signaling NaNs
2551	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
2552	inline half operator+(half x, half y)
2553	{
2554	#ifdef HALF_ARITHMETIC_TYPE
2555		return half(detail::binary, detail::float2half<half::round_style>(detail::half2float<detail::internal_t>(x.data_)+detail::half2float<detail::internal_t>(y.data_)));
2556	#else
2557		int absx = x.data_ & 0x7FFF, absy = y.data_ & 0x7FFF;
2558		bool sub = ((x.data_^y.data_)&0x8000) != 0;
2559		if(absx >= 0x7C00 || absy >= 0x7C00)
2560			return half(detail::binary,	(absx>0x7C00 || absy>0x7C00) ? detail::signal(x.data_, y.data_) : (absy!=0x7C00) ? x.data_ :
2561										(sub && absx==0x7C00) ? detail::invalid() : y.data_);
2562		if(!absx)
2563			return absy ? y : half(detail::binary, (half::round_style==std::round_toward_neg_infinity) ? (x.data_|y.data_) : (x.data_&y.data_));
2564		if(!absy)
2565			return x;
2566		unsigned int sign = ((sub && absy>absx) ? y.data_ : x.data_) & 0x8000;
2567		if(absy > absx)
2568			std::swap(absx, absy);
2569		int exp = (absx>>10) + (absx<=0x3FF), d = exp - (absy>>10) - (absy<=0x3FF), mx = ((absx&0x3FF)|((absx>0x3FF)<<10)) << 3, my;
2570		if(d < 13)
2571		{
2572			my = ((absy&0x3FF)|((absy>0x3FF)<<10)) << 3;
2573			my = (my>>d) | ((my&((1<<d)-1))!=0);
2574		}
2575		else
2576			my = 1;
2577		if(sub)
2578		{
2579			if(!(mx-=my))
2580				return half(detail::binary, static_cast<unsigned>(half::round_style==std::round_toward_neg_infinity)<<15);
2581			for(; mx<0x2000 && exp>1; mx<<=1,--exp) ;
2582		}
2583		else
2584		{
2585			mx += my;
2586			int i = mx >> 14;
2587			if((exp+=i) > 30)
2588				return half(detail::binary, detail::overflow<half::round_style>(sign));
2589			mx = (mx>>i) | (mx&i);
2590		}
2591		return half(detail::binary, detail::rounded<half::round_style,false>(sign+((exp-1)<<10)+(mx>>3), (mx>>2)&1, (mx&0x3)!=0));
2592	#endif
2593	}
2594
2595	/// Subtraction.
2596	/// This operation is exact to rounding for all rounding modes.
2597	/// \param x left operand
2598	/// \param y right operand
2599	/// \return difference of half expressions
2600	/// \exception FE_INVALID if \a x and \a y are infinities with equal signs or signaling NaNs
2601	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
2602	inline half operator-(half x, half y)
2603	{
2604	#ifdef HALF_ARITHMETIC_TYPE
2605		return half(detail::binary, detail::float2half<half::round_style>(detail::half2float<detail::internal_t>(x.data_)-detail::half2float<detail::internal_t>(y.data_)));
2606	#else
2607		return x + -y;
2608	#endif
2609	}
2610
2611	/// Multiplication.
2612	/// This operation is exact to rounding for all rounding modes.
2613	/// \param x left operand
2614	/// \param y right operand
2615	/// \return product of half expressions
2616	/// \exception FE_INVALID if multiplying 0 with infinity or if \a x or \a y is signaling NaN
2617	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
2618	inline half operator*(half x, half y)
2619	{
2620	#ifdef HALF_ARITHMETIC_TYPE
2621		return half(detail::binary, detail::float2half<half::round_style>(detail::half2float<detail::internal_t>(x.data_)*detail::half2float<detail::internal_t>(y.data_)));
2622	#else
2623		int absx = x.data_ & 0x7FFF, absy = y.data_ & 0x7FFF, exp = -16;
2624		unsigned int sign = (x.data_^y.data_) & 0x8000;
2625		if(absx >= 0x7C00 || absy >= 0x7C00)
2626			return half(detail::binary,	(absx>0x7C00 || absy>0x7C00) ? detail::signal(x.data_, y.data_) :
2627										((absx==0x7C00 && !absy)||(absy==0x7C00 && !absx)) ? detail::invalid() : (sign|0x7C00));
2628		if(!absx || !absy)
2629			return half(detail::binary, sign);
2630		for(; absx<0x400; absx<<=1,--exp) ;
2631		for(; absy<0x400; absy<<=1,--exp) ;
2632		detail::uint32 m = static_cast<detail::uint32>((absx&0x3FF)|0x400) * static_cast<detail::uint32>((absy&0x3FF)|0x400);
2633		int i = m >> 21, s = m & i;
2634		exp += (absx>>10) + (absy>>10) + i;
2635		if(exp > 29)
2636			return half(detail::binary, detail::overflow<half::round_style>(sign));
2637		else if(exp < -11)
2638			return half(detail::binary, detail::underflow<half::round_style>(sign));
2639		return half(detail::binary, detail::fixed2half<half::round_style,20,false,false,false>(m>>i, exp, sign, s));
2640	#endif
2641	}
2642
2643	/// Division.
2644	/// This operation is exact to rounding for all rounding modes.
2645	/// \param x left operand
2646	/// \param y right operand
2647	/// \return quotient of half expressions
2648	/// \exception FE_INVALID if dividing 0s or infinities with each other or if \a x or \a y is signaling NaN
2649	/// \exception FE_DIVBYZERO if dividing finite value by 0
2650	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
2651	inline half operator/(half x, half y)
2652	{
2653	#ifdef HALF_ARITHMETIC_TYPE
2654		return half(detail::binary, detail::float2half<half::round_style>(detail::half2float<detail::internal_t>(x.data_)/detail::half2float<detail::internal_t>(y.data_)));
2655	#else
2656		int absx = x.data_ & 0x7FFF, absy = y.data_ & 0x7FFF, exp = 14;
2657		unsigned int sign = (x.data_^y.data_) & 0x8000;
2658		if(absx >= 0x7C00 || absy >= 0x7C00)
2659			return half(detail::binary,	(absx>0x7C00 || absy>0x7C00) ? detail::signal(x.data_, y.data_) :
2660										(absx==absy) ? detail::invalid() : (sign|((absx==0x7C00) ? 0x7C00 : 0)));
2661		if(!absx)
2662			return half(detail::binary, absy ? sign : detail::invalid());
2663		if(!absy)
2664			return half(detail::binary, detail::pole(sign));
2665		for(; absx<0x400; absx<<=1,--exp) ;
2666		for(; absy<0x400; absy<<=1,++exp) ;
2667		detail::uint32 mx = (absx&0x3FF) | 0x400, my = (absy&0x3FF) | 0x400;
2668		int i = mx < my;
2669		exp += (absx>>10) - (absy>>10) - i;
2670		if(exp > 29)
2671			return half(detail::binary, detail::overflow<half::round_style>(sign));
2672		else if(exp < -11)
2673			return half(detail::binary, detail::underflow<half::round_style>(sign));
2674		mx <<= 12 + i;
2675		my <<= 1;
2676		return half(detail::binary, detail::fixed2half<half::round_style,11,false,false,false>(mx/my, exp, sign, mx%my!=0));
2677	#endif
2678	}
2679
2680	/// \}
2681	/// \anchor streaming
2682	/// \name Input and output
2683	/// \{
2684
2685	/// Output operator.
2686	///	This uses the built-in functionality for streaming out floating-point numbers.
2687	/// \param out output stream to write into
2688	/// \param arg half expression to write
2689	/// \return reference to output stream
2690	template<typename charT,typename traits> std::basic_ostream<charT,traits>& operator<<(std::basic_ostream<charT,traits> &out, half arg)
2691	{
2692	#ifdef HALF_ARITHMETIC_TYPE
2693		return out << detail::half2float<detail::internal_t>(arg.data_);
2694	#else
2695		return out << detail::half2float<float>(arg.data_);
2696	#endif
2697	}
2698
2699	/// Input operator.
2700	///	This uses the built-in functionality for streaming in floating-point numbers, specifically double precision floating 
2701	/// point numbers (unless overridden with [HALF_ARITHMETIC_TYPE](\ref HALF_ARITHMETIC_TYPE)). So the input string is first 
2702	/// rounded to double precision using the underlying platform's current floating-point rounding mode before being rounded 
2703	/// to half-precision using the library's half-precision rounding mode.
2704	/// \param in input stream to read from
2705	/// \param arg half to read into
2706	/// \return reference to input stream
2707	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
2708	template<typename charT,typename traits> std::basic_istream<charT,traits>& operator>>(std::basic_istream<charT,traits> &in, half &arg)
2709	{
2710	#ifdef HALF_ARITHMETIC_TYPE
2711		detail::internal_t f;
2712	#else
2713		double f;
2714	#endif
2715		if(in >> f)
2716			arg.data_ = detail::float2half<half::round_style>(f);
2717		return in;
2718	}
2719
2720	/// \}
2721	/// \anchor basic
2722	/// \name Basic mathematical operations
2723	/// \{
2724
2725	/// Absolute value.
2726	/// **See also:** Documentation for [std::fabs](https://en.cppreference.com/w/cpp/numeric/math/fabs).
2727	/// \param arg operand
2728	/// \return absolute value of \a arg
2729	inline HALF_CONSTEXPR half fabs(half arg) { return half(detail::binary, arg.data_&0x7FFF); }
2730
2731	/// Absolute value.
2732	/// **See also:** Documentation for [std::abs](https://en.cppreference.com/w/cpp/numeric/math/fabs).
2733	/// \param arg operand
2734	/// \return absolute value of \a arg
2735	inline HALF_CONSTEXPR half abs(half arg) { return fabs(arg); }
2736
2737	/// Remainder of division.
2738	/// **See also:** Documentation for [std::fmod](https://en.cppreference.com/w/cpp/numeric/math/fmod).
2739	/// \param x first operand
2740	/// \param y second operand
2741	/// \return remainder of floating-point division.
2742	/// \exception FE_INVALID if \a x is infinite or \a y is 0 or if \a x or \a y is signaling NaN
2743	inline half fmod(half x, half y)
2744	{
2745		unsigned int absx = x.data_ & 0x7FFF, absy = y.data_ & 0x7FFF, sign = x.data_ & 0x8000;
2746		if(absx >= 0x7C00 || absy >= 0x7C00)
2747			return half(detail::binary,	(absx>0x7C00 || absy>0x7C00) ? detail::signal(x.data_, y.data_) :
2748										(absx==0x7C00) ? detail::invalid() : x.data_);
2749		if(!absy)
2750			return half(detail::binary, detail::invalid());
2751		if(!absx)
2752			return x;
2753		if(absx == absy)
2754			return half(detail::binary, sign);
2755		return half(detail::binary, sign|detail::mod<false,false>(absx, absy));
2756	}
2757
2758	/// Remainder of division.
2759	/// **See also:** Documentation for [std::remainder](https://en.cppreference.com/w/cpp/numeric/math/remainder).
2760	/// \param x first operand
2761	/// \param y second operand
2762	/// \return remainder of floating-point division.
2763	/// \exception FE_INVALID if \a x is infinite or \a y is 0 or if \a x or \a y is signaling NaN
2764	inline half remainder(half x, half y)
2765	{
2766		unsigned int absx = x.data_ & 0x7FFF, absy = y.data_ & 0x7FFF, sign = x.data_ & 0x8000;
2767		if(absx >= 0x7C00 || absy >= 0x7C00)
2768			return half(detail::binary,	(absx>0x7C00 || absy>0x7C00) ? detail::signal(x.data_, y.data_) :
2769										(absx==0x7C00) ? detail::invalid() : x.data_);
2770		if(!absy)
2771			return half(detail::binary, detail::invalid());
2772		if(absx == absy)
2773			return half(detail::binary, sign);
2774		return half(detail::binary, sign^detail::mod<false,true>(absx, absy));
2775	}
2776
2777	/// Remainder of division.
2778	/// **See also:** Documentation for [std::remquo](https://en.cppreference.com/w/cpp/numeric/math/remquo).
2779	/// \param x first operand
2780	/// \param y second operand
2781	/// \param quo address to store some bits of quotient at
2782	/// \return remainder of floating-point division.
2783	/// \exception FE_INVALID if \a x is infinite or \a y is 0 or if \a x or \a y is signaling NaN
2784	inline half remquo(half x, half y, int *quo)
2785	{
2786		unsigned int absx = x.data_ & 0x7FFF, absy = y.data_ & 0x7FFF, value = x.data_ & 0x8000;
2787		if(absx >= 0x7C00 || absy >= 0x7C00)
2788			return half(detail::binary,	(absx>0x7C00 || absy>0x7C00) ? detail::signal(x.data_, y.data_) :
2789										(absx==0x7C00) ? detail::invalid() : (*quo = 0, x.data_));
2790		if(!absy)
2791			return half(detail::binary, detail::invalid());
2792		bool qsign = ((value^y.data_)&0x8000) != 0;
2793		int q = 1;
2794		if(absx != absy)
2795			value ^= detail::mod<true, true>(absx, absy, &q);
2796		return *quo = qsign ? -q : q, half(detail::binary, value);
2797	}
2798
2799	/// Fused multiply add.
2800	/// This function is exact to rounding for all rounding modes.
2801	///
2802	/// **See also:** Documentation for [std::fma](https://en.cppreference.com/w/cpp/numeric/math/fma).
2803	/// \param x first operand
2804	/// \param y second operand
2805	/// \param z third operand
2806	/// \return ( \a x * \a y ) + \a z rounded as one operation.
2807	/// \exception FE_INVALID according to operator*() and operator+() unless any argument is a quiet NaN and no argument is a signaling NaN
2808	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding the final addition
2809	inline half fma(half x, half y, half z)
2810	{
2811	#ifdef HALF_ARITHMETIC_TYPE
2812		detail::internal_t fx = detail::half2float<detail::internal_t>(x.data_), fy = detail::half2float<detail::internal_t>(y.data_), fz = detail::half2float<detail::internal_t>(z.data_);
2813		#if HALF_ENABLE_CPP11_CMATH && FP_FAST_FMA
2814			return half(detail::binary, detail::float2half<half::round_style>(std::fma(fx, fy, fz)));
2815		#else
2816			return half(detail::binary, detail::float2half<half::round_style>(fx*fy+fz));
2817		#endif
2818	#else
2819		int absx = x.data_ & 0x7FFF, absy = y.data_ & 0x7FFF, absz = z.data_ & 0x7FFF, exp = -15;
2820		unsigned int sign = (x.data_^y.data_) & 0x8000;
2821		bool sub = ((sign^z.data_)&0x8000) != 0;
2822		if(absx >= 0x7C00 || absy >= 0x7C00 || absz >= 0x7C00)
2823			return	(absx>0x7C00 || absy>0x7C00 || absz>0x7C00) ? half(detail::binary, detail::signal(x.data_, y.data_, z.data_)) :
2824					(absx==0x7C00) ? half(detail::binary, (!absy || (sub && absz==0x7C00)) ? detail::invalid() : (sign|0x7C00)) :
2825					(absy==0x7C00) ? half(detail::binary, (!absx || (sub && absz==0x7C00)) ? detail::invalid() : (sign|0x7C00)) : z;
2826		if(!absx || !absy)
2827			return absz ? z : half(detail::binary, (half::round_style==std::round_toward_neg_infinity) ? (z.data_|sign) : (z.data_&sign));
2828		for(; absx<0x400; absx<<=1,--exp) ;
2829		for(; absy<0x400; absy<<=1,--exp) ;
2830		detail::uint32 m = static_cast<detail::uint32>((absx&0x3FF)|0x400) * static_cast<detail::uint32>((absy&0x3FF)|0x400);
2831		int i = m >> 21;
2832		exp += (absx>>10) + (absy>>10) + i;
2833		m <<= 3 - i;
2834		if(absz)
2835		{
2836			int expz = 0;
2837			for(; absz<0x400; absz<<=1,--expz) ;
2838			expz += absz >> 10;
2839			detail::uint32 mz = static_cast<detail::uint32>((absz&0x3FF)|0x400) << 13;
2840			if(expz > exp || (expz == exp && mz > m))
2841			{
2842				std::swap(m, mz);
2843				std::swap(exp, expz);
2844				if(sub)
2845					sign = z.data_ & 0x8000;
2846			}
2847			int d = exp - expz;
2848			mz = (d<23) ? ((mz>>d)|((mz&((static_cast<detail::uint32>(1)<<d)-1))!=0)) : 1;
2849			if(sub)
2850			{
2851				m = m - mz;
2852				if(!m)
2853					return half(detail::binary, static_cast<unsigned>(half::round_style==std::round_toward_neg_infinity)<<15);
2854				for(; m<0x800000; m<<=1,--exp) ;
2855			}
2856			else
2857			{
2858				m += mz;
2859				i = m >> 24;
2860				m = (m>>i) | (m&i);
2861				exp += i;
2862			}
2863		}
2864		if(exp > 30)
2865			return half(detail::binary, detail::overflow<half::round_style>(sign));
2866		else if(exp < -10)
2867			return half(detail::binary, detail::underflow<half::round_style>(sign));
2868		return half(detail::binary, detail::fixed2half<half::round_style,23,false,false,false>(m, exp-1, sign));
2869	#endif
2870	}
2871
2872	/// Maximum of half expressions.
2873	/// **See also:** Documentation for [std::fmax](https://en.cppreference.com/w/cpp/numeric/math/fmax).
2874	/// \param x first operand
2875	/// \param y second operand
2876	/// \return maximum of operands, ignoring quiet NaNs
2877	/// \exception FE_INVALID if \a x or \a y is signaling NaN
2878	inline HALF_CONSTEXPR_NOERR half fmax(half x, half y)
2879	{
2880		return half(detail::binary, (!isnan(y) && (isnan(x) || (x.data_^(0x8000|(0x8000-(x.data_>>15)))) < 
2881			(y.data_^(0x8000|(0x8000-(y.data_>>15)))))) ? detail::select(y.data_, x.data_) : detail::select(x.data_, y.data_));
2882	}
2883
2884	/// Minimum of half expressions.
2885	/// **See also:** Documentation for [std::fmin](https://en.cppreference.com/w/cpp/numeric/math/fmin).
2886	/// \param x first operand
2887	/// \param y second operand
2888	/// \return minimum of operands, ignoring quiet NaNs
2889	/// \exception FE_INVALID if \a x or \a y is signaling NaN
2890	inline HALF_CONSTEXPR_NOERR half fmin(half x, half y)
2891	{
2892		return half(detail::binary, (!isnan(y) && (isnan(x) || (x.data_^(0x8000|(0x8000-(x.data_>>15)))) >
2893			(y.data_^(0x8000|(0x8000-(y.data_>>15)))))) ? detail::select(y.data_, x.data_) : detail::select(x.data_, y.data_));
2894	}
2895
2896	/// Positive difference.
2897	/// This function is exact to rounding for all rounding modes.
2898	///
2899	/// **See also:** Documentation for [std::fdim](https://en.cppreference.com/w/cpp/numeric/math/fdim).
2900	/// \param x first operand
2901	/// \param y second operand
2902	/// \return \a x - \a y or 0 if difference negative
2903	/// \exception FE_... according to operator-(half,half)
2904	inline half fdim(half x, half y)
2905	{
2906		if(isnan(x) || isnan(y))
2907			return half(detail::binary, detail::signal(x.data_, y.data_));
2908		return (x.data_^(0x8000|(0x8000-(x.data_>>15)))) <= (y.data_^(0x8000|(0x8000-(y.data_>>15)))) ? half(detail::binary, 0) : (x-y);
2909	}
2910
2911	/// Get NaN value.
2912	/// **See also:** Documentation for [std::nan](https://en.cppreference.com/w/cpp/numeric/math/nan).
2913	/// \param arg string code
2914	/// \return quiet NaN
2915	inline half nanh(const char *arg)
2916	{
2917		unsigned int value = 0x7FFF;
2918		while(*arg)
2919			value ^= static_cast<unsigned>(*arg++) & 0xFF;
2920		return half(detail::binary, value);
2921	}
2922
2923	/// \}
2924	/// \anchor exponential
2925	/// \name Exponential functions
2926	/// \{
2927
2928	/// Exponential function.
2929	/// This function is exact to rounding for all rounding modes.
2930	///
2931	/// **See also:** Documentation for [std::exp](https://en.cppreference.com/w/cpp/numeric/math/exp).
2932	/// \param arg function argument
2933	/// \return e raised to \a arg
2934	/// \exception FE_INVALID for signaling NaN
2935	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
2936	inline half exp(half arg)
2937	{
2938	#ifdef HALF_ARITHMETIC_TYPE
2939		return half(detail::binary, detail::float2half<half::round_style>(std::exp(detail::half2float<detail::internal_t>(arg.data_))));
2940	#else
2941		int abs = arg.data_ & 0x7FFF, e = (abs>>10) + (abs<=0x3FF), exp;
2942		if(!abs)
2943			return half(detail::binary, 0x3C00);
2944		if(abs >= 0x7C00)
2945			return half(detail::binary, (abs==0x7C00) ? (0x7C00&((arg.data_>>15)-1U)) : detail::signal(arg.data_));
2946		if(abs >= 0x4C80)
2947			return half(detail::binary, (arg.data_&0x8000) ? detail::underflow<half::round_style>() : detail::overflow<half::round_style>());
2948		detail::uint32 m = detail::multiply64(static_cast<detail::uint32>((abs&0x3FF)+((abs>0x3FF)<<10))<<21, 0xB8AA3B29);
2949		if(e < 14)
2950		{
2951			exp = 0;
2952			m >>= 14 - e;
2953		}
2954		else
2955		{
2956			exp = m >> (45-e);
2957			m = (m<<(e-14)) & 0x7FFFFFFF;
2958		}
2959		return half(detail::binary, detail::exp2_post<half::round_style>(m, exp, (arg.data_&0x8000)!=0, 0, 26));
2960	#endif
2961	}
2962
2963	/// Binary exponential.
2964	/// This function is exact to rounding for all rounding modes.
2965	///
2966	/// **See also:** Documentation for [std::exp2](https://en.cppreference.com/w/cpp/numeric/math/exp2).
2967	/// \param arg function argument
2968	/// \return 2 raised to \a arg
2969	/// \exception FE_INVALID for signaling NaN
2970	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
2971	inline half exp2(half arg)
2972	{
2973	#if defined(HALF_ARITHMETIC_TYPE) && HALF_ENABLE_CPP11_CMATH
2974		return half(detail::binary, detail::float2half<half::round_style>(std::exp2(detail::half2float<detail::internal_t>(arg.data_))));
2975	#else
2976		int abs = arg.data_ & 0x7FFF, e = (abs>>10) + (abs<=0x3FF), exp = (abs&0x3FF) + ((abs>0x3FF)<<10);
2977		if(!abs)
2978			return half(detail::binary, 0x3C00);
2979		if(abs >= 0x7C00)
2980			return half(detail::binary, (abs==0x7C00) ? (0x7C00&((arg.data_>>15)-1U)) : detail::signal(arg.data_));
2981		if(abs >= 0x4E40)
2982			return half(detail::binary, (arg.data_&0x8000) ? detail::underflow<half::round_style>() : detail::overflow<half::round_style>());
2983		return half(detail::binary, detail::exp2_post<half::round_style>(
2984			(static_cast<detail::uint32>(exp)<<(6+e))&0x7FFFFFFF, exp>>(25-e), (arg.data_&0x8000)!=0, 0, 28));
2985	#endif
2986	}
2987
2988	/// Exponential minus one.
2989	/// This function may be 1 ULP off the correctly rounded exact result in <0.05% of inputs for `std::round_to_nearest` 
2990	/// and in <1% of inputs for any other rounding mode.
2991	///
2992	/// **See also:** Documentation for [std::expm1](https://en.cppreference.com/w/cpp/numeric/math/expm1).
2993	/// \param arg function argument
2994	/// \return e raised to \a arg and subtracted by 1
2995	/// \exception FE_INVALID for signaling NaN
2996	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
2997	inline half expm1(half arg)
2998	{
2999	#if defined(HALF_ARITHMETIC_TYPE) && HALF_ENABLE_CPP11_CMATH
3000		return half(detail::binary, detail::float2half<half::round_style>(std::expm1(detail::half2float<detail::internal_t>(arg.data_))));
3001	#else
3002		unsigned int abs = arg.data_ & 0x7FFF, sign = arg.data_ & 0x8000, e = (abs>>10) + (abs<=0x3FF), exp;
3003		if(!abs)
3004			return arg;
3005		if(abs >= 0x7C00)
3006			return half(detail::binary, (abs==0x7C00) ? (0x7C00+(sign>>1)) : detail::signal(arg.data_));
3007		if(abs >= 0x4A00)
3008			return half(detail::binary, (arg.data_&0x8000) ? detail::rounded<half::round_style,true>(0xBBFF, 1, 1) : detail::overflow<half::round_style>());
3009		detail::uint32 m = detail::multiply64(static_cast<detail::uint32>((abs&0x3FF)+((abs>0x3FF)<<10))<<21, 0xB8AA3B29);
3010		if(e < 14)
3011		{
3012			exp = 0;
3013			m >>= 14 - e;
3014		}
3015		else
3016		{
3017			exp = m >> (45-e);
3018			m = (m<<(e-14)) & 0x7FFFFFFF;
3019		}
3020		m = detail::exp2(m);
3021		if(sign)
3022		{
3023			int s = 0;
3024			if(m > 0x80000000)
3025			{
3026				++exp;
3027				m = detail::divide64(0x80000000, m, s);
3028			}
3029			m = 0x80000000 - ((m>>exp)|((m&((static_cast<detail::uint32>(1)<<exp)-1))!=0)|s);
3030			exp = 0;
3031		}
3032		else
3033			m -= (exp<31) ? (0x80000000>>exp) : 1;
3034		for(exp+=14; m<0x80000000 && exp; m<<=1,--exp) ;
3035		if(exp > 29)
3036			return half(detail::binary, detail::overflow<half::round_style>());
3037		return half(detail::binary, detail::rounded<half::round_style,true>(sign+(exp<<10)+(m>>21), (m>>20)&1, (m&0xFFFFF)!=0));
3038	#endif
3039	}
3040
3041	/// Natural logarithm.
3042	/// This function is exact to rounding for all rounding modes.
3043	///
3044	/// **See also:** Documentation for [std::log](https://en.cppreference.com/w/cpp/numeric/math/log).
3045	/// \param arg function argument
3046	/// \return logarithm of \a arg to base e
3047	/// \exception FE_INVALID for signaling NaN or negative argument
3048	/// \exception FE_DIVBYZERO for 0
3049	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
3050	inline half log(half arg)
3051	{
3052	#ifdef HALF_ARITHMETIC_TYPE
3053		return half(detail::binary, detail::float2half<half::round_style>(std::log(detail::half2float<detail::internal_t>(arg.data_))));
3054	#else
3055		int abs = arg.data_ & 0x7FFF, exp = -15;
3056		if(!abs)
3057			return half(detail::binary, detail::pole(0x8000));
3058		if(arg.data_ & 0x8000)
3059			return half(detail::binary, (arg.data_<=0xFC00) ? detail::invalid() : detail::signal(arg.data_));
3060		if(abs >= 0x7C00)
3061			return (abs==0x7C00) ? arg : half(detail::binary, detail::signal(arg.data_));
3062		for(; abs<0x400; abs<<=1,--exp) ;
3063		exp += abs >> 10;
3064		return half(detail::binary, detail::log2_post<half::round_style,0xB8AA3B2A>(
3065			detail::log2(static_cast<detail::uint32>((abs&0x3FF)|0x400)<<20, 27)+8, exp, 17));
3066	#endif
3067	}
3068
3069	/// Common logarithm.
3070	/// This function is exact to rounding for all rounding modes.
3071	///
3072	/// **See also:** Documentation for [std::log10](https://en.cppreference.com/w/cpp/numeric/math/log10).
3073	/// \param arg function argument
3074	/// \return logarithm of \a arg to base 10
3075	/// \exception FE_INVALID for signaling NaN or negative argument
3076	/// \exception FE_DIVBYZERO for 0
3077	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
3078	inline half log10(half arg)
3079	{
3080	#ifdef HALF_ARITHMETIC_TYPE
3081		return half(detail::binary, detail::float2half<half::round_style>(std::log10(detail::half2float<detail::internal_t>(arg.data_))));
3082	#else
3083		int abs = arg.data_ & 0x7FFF, exp = -15;
3084		if(!abs)
3085			return half(detail::binary, detail::pole(0x8000));
3086		if(arg.data_ & 0x8000)
3087			return half(detail::binary, (arg.data_<=0xFC00) ? detail::invalid() : detail::signal(arg.data_));
3088		if(abs >= 0x7C00)
3089			return (abs==0x7C00) ? arg : half(detail::binary, detail::signal(arg.data_));
3090		switch(abs)
3091		{
3092			case 0x4900: return half(detail::binary, 0x3C00);
3093			case 0x5640: return half(detail::binary, 0x4000);
3094			case 0x63D0: return half(detail::binary, 0x4200);
3095			case 0x70E2: return half(detail::binary, 0x4400);
3096		}
3097		for(; abs<0x400; abs<<=1,--exp) ;
3098		exp += abs >> 10;
3099		return half(detail::binary, detail::log2_post<half::round_style,0xD49A784C>(
3100			detail::log2(static_cast<detail::uint32>((abs&0x3FF)|0x400)<<20, 27)+8, exp, 16));
3101	#endif
3102	}
3103
3104	/// Binary logarithm.
3105	/// This function is exact to rounding for all rounding modes.
3106	///
3107	/// **See also:** Documentation for [std::log2](https://en.cppreference.com/w/cpp/numeric/math/log2).
3108	/// \param arg function argument
3109	/// \return logarithm of \a arg to base 2
3110	/// \exception FE_INVALID for signaling NaN or negative argument
3111	/// \exception FE_DIVBYZERO for 0
3112	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
3113	inline half log2(half arg)
3114	{
3115	#if defined(HALF_ARITHMETIC_TYPE) && HALF_ENABLE_CPP11_CMATH
3116		return half(detail::binary, detail::float2half<half::round_style>(std::log2(detail::half2float<detail::internal_t>(arg.data_))));
3117	#else
3118		int abs = arg.data_ & 0x7FFF, exp = -15, s = 0;
3119		if(!abs)
3120			return half(detail::binary, detail::pole(0x8000));
3121		if(arg.data_ & 0x8000)
3122			return half(detail::binary, (arg.data_<=0xFC00) ? detail::invalid() : detail::signal(arg.data_));
3123		if(abs >= 0x7C00)
3124			return (abs==0x7C00) ? arg : half(detail::binary, detail::signal(arg.data_));
3125		if(abs == 0x3C00)
3126			return half(detail::binary, 0);
3127		for(; abs<0x400; abs<<=1,--exp) ;
3128		exp += (abs>>10);
3129		if(!(abs&0x3FF))
3130		{
3131			unsigned int value = static_cast<unsigned>(exp<0) << 15, m = std::abs(exp) << 6;
3132			for(exp=18; m<0x400; m<<=1,--exp) ;
3133			return half(detail::binary, value+(exp<<10)+m);
3134		}
3135		detail::uint32 ilog = exp, sign = detail::sign_mask(ilog), m = 
3136			(((ilog<<27)+(detail::log2(static_cast<detail::uint32>((abs&0x3FF)|0x400)<<20, 28)>>4))^sign) - sign;
3137		if(!m)
3138			return half(detail::binary, 0);
3139		for(exp=14; m<0x8000000 && exp; m<<=1,--exp) ;
3140		for(; m>0xFFFFFFF; m>>=1,++exp)
3141			s |= m & 1;
3142		return half(detail::binary, detail::fixed2half<half::round_style,27,false,false,true>(m, exp, sign&0x8000, s));
3143	#endif
3144	}
3145
3146	/// Natural logarithm plus one.
3147	/// This function may be 1 ULP off the correctly rounded exact result in <0.05% of inputs for `std::round_to_nearest` 
3148	/// and in ~1% of inputs for any other rounding mode.
3149	///
3150	/// **See also:** Documentation for [std::log1p](https://en.cppreference.com/w/cpp/numeric/math/log1p).
3151	/// \param arg function argument
3152	/// \return logarithm of \a arg plus 1 to base e
3153	/// \exception FE_INVALID for signaling NaN or argument <-1
3154	/// \exception FE_DIVBYZERO for -1
3155	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
3156	inline half log1p(half arg)
3157	{
3158	#if defined(HALF_ARITHMETIC_TYPE) && HALF_ENABLE_CPP11_CMATH
3159		return half(detail::binary, detail::float2half<half::round_style>(std::log1p(detail::half2float<detail::internal_t>(arg.data_))));
3160	#else
3161		if(arg.data_ >= 0xBC00)
3162			return half(detail::binary, (arg.data_==0xBC00) ? detail::pole(0x8000) : (arg.data_<=0xFC00) ? detail::invalid() : detail::signal(arg.data_));
3163		int abs = arg.data_ & 0x7FFF, exp = -15;
3164		if(!abs || abs >= 0x7C00)
3165			return (abs>0x7C00) ? half(detail::binary, detail::signal(arg.data_)) : arg;
3166		for(; abs<0x400; abs<<=1,--exp) ;
3167		exp += abs >> 10;
3168		detail::uint32 m = static_cast<detail::uint32>((abs&0x3FF)|0x400) << 20;
3169		if(arg.data_ & 0x8000)
3170		{
3171			m = 0x40000000 - (m>>-exp);
3172			for(exp=0; m<0x40000000; m<<=1,--exp) ;
3173		}
3174		else
3175		{
3176			if(exp < 0)
3177			{
3178				m = 0x40000000 + (m>>-exp);
3179				exp = 0;
3180			}
3181			else
3182			{
3183				m += 0x40000000 >> exp;
3184				int i = m >> 31;
3185				m >>= i;
3186				exp += i;
3187			}
3188		}
3189		return half(detail::binary, detail::log2_post<half::round_style,0xB8AA3B2A>(detail::log2(m), exp, 17));
3190	#endif
3191	}
3192
3193	/// \}
3194	/// \anchor power
3195	/// \name Power functions
3196	/// \{
3197
3198	/// Square root.
3199	/// This function is exact to rounding for all rounding modes.
3200	///
3201	/// **See also:** Documentation for [std::sqrt](https://en.cppreference.com/w/cpp/numeric/math/sqrt).
3202	/// \param arg function argument
3203	/// \return square root of \a arg
3204	/// \exception FE_INVALID for signaling NaN and negative arguments
3205	/// \exception FE_INEXACT according to rounding
3206	inline half sqrt(half arg)
3207	{
3208	#ifdef HALF_ARITHMETIC_TYPE
3209		return half(detail::binary, detail::float2half<half::round_style>(std::sqrt(detail::half2float<detail::internal_t>(arg.data_))));
3210	#else
3211		int abs = arg.data_ & 0x7FFF, exp = 15;
3212		if(!abs || arg.data_ >= 0x7C00)
3213			return half(detail::binary, (abs>0x7C00) ? detail::signal(arg.data_) : (arg.data_>0x8000) ? detail::invalid() : arg.data_);
3214		for(; abs<0x400; abs<<=1,--exp) ;
3215		detail::uint32 r = static_cast<detail::uint32>((abs&0x3FF)|0x400) << 10, m = detail::sqrt<20>(r, exp+=abs>>10);
3216		return half(detail::binary, detail::rounded<half::round_style,false>((exp<<10)+(m&0x3FF), r>m, r!=0));
3217	#endif
3218	}
3219
3220	/// Inverse square root.
3221	/// This function is exact to rounding for all rounding modes and thus generally more accurate than directly computing 
3222	/// 1 / sqrt(\a arg) in half-precision, in addition to also being faster.
3223	/// \param arg function argument
3224	/// \return reciprocal of square root of \a arg
3225	/// \exception FE_INVALID for signaling NaN and negative arguments
3226	/// \exception FE_INEXACT according to rounding
3227	inline half rsqrt(half arg)
3228	{
3229	#ifdef HALF_ARITHMETIC_TYPE
3230		return half(detail::binary, detail::float2half<half::round_style>(detail::internal_t(1)/std::sqrt(detail::half2float<detail::internal_t>(arg.data_))));
3231	#else
3232		unsigned int abs = arg.data_ & 0x7FFF, bias = 0x4000;
3233		if(!abs || arg.data_ >= 0x7C00)
3234			return half(detail::binary,	(abs>0x7C00) ? detail::signal(arg.data_) : (arg.data_>0x8000) ?
3235										detail::invalid() : !abs ? detail::pole(arg.data_&0x8000) : 0);
3236		for(; abs<0x400; abs<<=1,bias-=0x400) ;
3237		unsigned int frac = (abs+=bias) & 0x7FF;
3238		if(frac == 0x400)
3239			return half(detail::binary, 0x7A00-(abs>>1));
3240		if((half::round_style == std::round_to_nearest && (frac == 0x3FE || frac == 0x76C)) ||
3241		   (half::round_style != std::round_to_nearest && (frac == 0x15A || frac == 0x3FC || frac == 0x401 || frac == 0x402 || frac == 0x67B)))
3242			return pow(arg, half(detail::binary, 0xB800));
3243		detail::uint32 f = 0x17376 - abs, mx = (abs&0x3FF) | 0x400, my = ((f>>1)&0x3FF) | 0x400, mz = my * my;
3244		int expy = (f>>11) - 31, expx = 32 - (abs>>10), i = mz >> 21;
3245		for(mz=0x60000000-(((mz>>i)*mx)>>(expx-2*expy-i)); mz<0x40000000; mz<<=1,--expy) ;
3246		i = (my*=mz>>10) >> 31;
3247		expy += i;
3248		my = (my>>(20+i)) + 1;
3249		i = (mz=my*my) >> 21;
3250		for(mz=0x60000000-(((mz>>i)*mx)>>(expx-2*expy-i)); mz<0x40000000; mz<<=1,--expy) ;
3251		i = (my*=(mz>>10)+1) >> 31;
3252		return half(detail::binary, detail::fixed2half<half::round_style,30,false,false,true>(my>>i, expy+i+14));
3253	#endif
3254	}
3255
3256	/// Cubic root.
3257	/// This function is exact to rounding for all rounding modes.
3258	///
3259	/// **See also:** Documentation for [std::cbrt](https://en.cppreference.com/w/cpp/numeric/math/cbrt).
3260	/// \param arg function argument
3261	/// \return cubic root of \a arg
3262	/// \exception FE_INVALID for signaling NaN
3263	/// \exception FE_INEXACT according to rounding
3264	inline half cbrt(half arg)
3265	{
3266	#if defined(HALF_ARITHMETIC_TYPE) && HALF_ENABLE_CPP11_CMATH
3267		return half(detail::binary, detail::float2half<half::round_style>(std::cbrt(detail::half2float<detail::internal_t>(arg.data_))));
3268	#else
3269		int abs = arg.data_ & 0x7FFF, exp = -15;
3270		if(!abs || abs == 0x3C00 || abs >= 0x7C00)
3271			return (abs>0x7C00) ? half(detail::binary, detail::signal(arg.data_)) : arg;
3272		for(; abs<0x400; abs<<=1, --exp);
3273		detail::uint32 ilog = exp + (abs>>10), sign = detail::sign_mask(ilog), f, m = 
3274			(((ilog<<27)+(detail::log2(static_cast<detail::uint32>((abs&0x3FF)|0x400)<<20, 24)>>4))^sign) - sign;
3275		for(exp=2; m<0x80000000; m<<=1,--exp) ;
3276		m = detail::multiply64(m, 0xAAAAAAAB);
3277		int i = m >> 31, s;
3278		exp += i;
3279		m <<= 1 - i;
3280		if(exp < 0)
3281		{
3282			f = m >> -exp;
3283			exp = 0;
3284		}
3285		else
3286		{
3287			f = (m<<exp) & 0x7FFFFFFF;
3288			exp = m >> (31-exp);
3289		}
3290		m = detail::exp2(f, (half::round_style==std::round_to_nearest) ? 29 : 26);
3291		if(sign)
3292		{
3293			if(m > 0x80000000)
3294			{
3295				m = detail::divide64(0x80000000, m, s);
3296				++exp;
3297			}
3298			exp = -exp;
3299		}
3300		return half(detail::binary, (half::round_style==std::round_to_nearest) ?
3301			detail::fixed2half<half::round_style,31,false,false,false>(m, exp+14, arg.data_&0x8000) :
3302			detail::fixed2half<half::round_style,23,false,false,false>((m+0x80)>>8, exp+14, arg.data_&0x8000));
3303	#endif
3304	}
3305
3306	/// Hypotenuse function.
3307	/// This function is exact to rounding for all rounding modes.
3308	///
3309	/// **See also:** Documentation for [std::hypot](https://en.cppreference.com/w/cpp/numeric/math/hypot).
3310	/// \param x first argument
3311	/// \param y second argument
3312	/// \return square root of sum of squares without internal over- or underflows
3313	/// \exception FE_INVALID if \a x or \a y is signaling NaN
3314	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding of the final square root
3315	inline half hypot(half x, half y)
3316	{
3317	#ifdef HALF_ARITHMETIC_TYPE
3318		detail::internal_t fx = detail::half2float<detail::internal_t>(x.data_), fy = detail::half2float<detail::internal_t>(y.data_);
3319		#if HALF_ENABLE_CPP11_CMATH
3320			return half(detail::binary, detail::float2half<half::round_style>(std::hypot(fx, fy)));
3321		#else
3322			return half(detail::binary, detail::float2half<half::round_style>(std::sqrt(fx*fx+fy*fy)));
3323		#endif
3324	#else
3325		int absx = x.data_ & 0x7FFF, absy = y.data_ & 0x7FFF, expx = 0, expy = 0;
3326		if(absx >= 0x7C00 || absy >= 0x7C00)
3327			return half(detail::binary,	(absx==0x7C00) ? detail::select(0x7C00, y.data_) :
3328				(absy==0x7C00) ? detail::select(0x7C00, x.data_) : detail::signal(x.data_, y.data_));
3329		if(!absx)
3330			return half(detail::binary, absy ? detail::check_underflow(absy) : 0);
3331		if(!absy)
3332			return half(detail::binary, detail::check_underflow(absx));
3333		if(absy > absx)
3334			std::swap(absx, absy);
3335		for(; absx<0x400; absx<<=1,--expx) ;
3336		for(; absy<0x400; absy<<=1,--expy) ;
3337		detail::uint32 mx = (absx&0x3FF) | 0x400, my = (absy&0x3FF) | 0x400;
3338		mx *= mx;
3339		my *= my;
3340		int ix = mx >> 21, iy = my >> 21;
3341		expx = 2*(expx+(absx>>10)) - 15 + ix;
3342		expy = 2*(expy+(absy>>10)) - 15 + iy;
3343		mx <<= 10 - ix;
3344		my <<= 10 - iy;
3345		int d = expx - expy;
3346		my = (d<30) ? ((my>>d)|((my&((static_cast<detail::uint32>(1)<<d)-1))!=0)) : 1;
3347		return half(detail::binary, detail::hypot_post<half::round_style>(mx+my, expx));
3348	#endif
3349	}
3350
3351	/// Hypotenuse function.
3352	/// This function is exact to rounding for all rounding modes.
3353	///
3354	/// **See also:** Documentation for [std::hypot](https://en.cppreference.com/w/cpp/numeric/math/hypot).
3355	/// \param x first argument
3356	/// \param y second argument
3357	/// \param z third argument
3358	/// \return square root of sum of squares without internal over- or underflows
3359	/// \exception FE_INVALID if \a x, \a y or \a z is signaling NaN
3360	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding of the final square root
3361	inline half hypot(half x, half y, half z)
3362	{
3363	#ifdef HALF_ARITHMETIC_TYPE
3364		detail::internal_t fx = detail::half2float<detail::internal_t>(x.data_), fy = detail::half2float<detail::internal_t>(y.data_), fz = detail::half2float<detail::internal_t>(z.data_);
3365		return half(detail::binary, detail::float2half<half::round_style>(std::sqrt(fx*fx+fy*fy+fz*fz)));
3366	#else
3367		int absx = x.data_ & 0x7FFF, absy = y.data_ & 0x7FFF, absz = z.data_ & 0x7FFF, expx = 0, expy = 0, expz = 0;
3368		if(!absx)
3369			return hypot(y, z);
3370		if(!absy)
3371			return hypot(x, z);
3372		if(!absz)
3373			return hypot(x, y);
3374		if(absx >= 0x7C00 || absy >= 0x7C00 || absz >= 0x7C00)
3375			return half(detail::binary,	(absx==0x7C00) ? detail::select(0x7C00, detail::select(y.data_, z.data_)) :
3376										(absy==0x7C00) ? detail::select(0x7C00, detail::select(x.data_, z.data_)) :
3377										(absz==0x7C00) ? detail::select(0x7C00, detail::select(x.data_, y.data_)) :
3378										detail::signal(x.data_, y.data_, z.data_));
3379		if(absz > absy)
3380			std::swap(absy, absz);
3381		if(absy > absx)
3382			std::swap(absx, absy);
3383		if(absz > absy)
3384			std::swap(absy, absz);
3385		for(; absx<0x400; absx<<=1,--expx) ;
3386		for(; absy<0x400; absy<<=1,--expy) ;
3387		for(; absz<0x400; absz<<=1,--expz) ;
3388		detail::uint32 mx = (absx&0x3FF) | 0x400, my = (absy&0x3FF) | 0x400, mz = (absz&0x3FF) | 0x400;
3389		mx *= mx;
3390		my *= my;
3391		mz *= mz;
3392		int ix = mx >> 21, iy = my >> 21, iz = mz >> 21;
3393		expx = 2*(expx+(absx>>10)) - 15 + ix;
3394		expy = 2*(expy+(absy>>10)) - 15 + iy;
3395		expz = 2*(expz+(absz>>10)) - 15 + iz;
3396		mx <<= 10 - ix;
3397		my <<= 10 - iy;
3398		mz <<= 10 - iz;
3399		int d = expy - expz;
3400		mz = (d<30) ? ((mz>>d)|((mz&((static_cast<detail::uint32>(1)<<d)-1))!=0)) : 1;
3401		my += mz;
3402		if(my & 0x80000000)
3403		{
3404			my = (my>>1) | (my&1);
3405			if(++expy > expx)
3406			{
3407				std::swap(mx, my);
3408				std::swap(expx, expy);
3409			}
3410		}
3411		d = expx - expy;
3412		my = (d<30) ? ((my>>d)|((my&((static_cast<detail::uint32>(1)<<d)-1))!=0)) : 1;
3413		return half(detail::binary, detail::hypot_post<half::round_style>(mx+my, expx));
3414	#endif
3415	}
3416
3417	/// Power function.
3418	/// This function may be 1 ULP off the correctly rounded exact result for any rounding mode in ~0.00025% of inputs.
3419	///
3420	/// **See also:** Documentation for [std::pow](https://en.cppreference.com/w/cpp/numeric/math/pow).
3421	/// \param x base
3422	/// \param y exponent
3423	/// \return \a x raised to \a y
3424	/// \exception FE_INVALID if \a x or \a y is signaling NaN or if \a x is finite an negative and \a y is finite and not integral
3425	/// \exception FE_DIVBYZERO if \a x is 0 and \a y is negative
3426	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
3427	inline half pow(half x, half y)
3428	{
3429	#ifdef HALF_ARITHMETIC_TYPE
3430		return half(detail::binary, detail::float2half<half::round_style>(std::pow(detail::half2float<detail::internal_t>(x.data_), detail::half2float<detail::internal_t>(y.data_))));
3431	#else
3432		int absx = x.data_ & 0x7FFF, absy = y.data_ & 0x7FFF, exp = -15;
3433		if(!absy || x.data_ == 0x3C00)
3434			return half(detail::binary, detail::select(0x3C00, (x.data_==0x3C00) ? y.data_ : x.data_));
3435		bool is_int = absy >= 0x6400 || (absy>=0x3C00 && !(absy&((1<<(25-(absy>>10)))-1)));
3436		unsigned int sign = x.data_ & (static_cast<unsigned>((absy<0x6800)&&is_int&&((absy>>(25-(absy>>10)))&1))<<15);
3437		if(absx >= 0x7C00 || absy >= 0x7C00)
3438			return half(detail::binary,	(absx>0x7C00 || absy>0x7C00) ? detail::signal(x.data_, y.data_) :
3439										(absy==0x7C00) ? ((absx==0x3C00) ? 0x3C00 : (!absx && y.data_==0xFC00) ? detail::pole() :
3440										(0x7C00&-((y.data_>>15)^(absx>0x3C00)))) : (sign|(0x7C00&((y.data_>>15)-1U))));
3441		if(!absx)
3442			return half(detail::binary, (y.data_&0x8000) ? detail::pole(sign) : sign);
3443		if((x.data_&0x8000) && !is_int)
3444			return half(detail::binary, detail::invalid());
3445		if(x.data_ == 0xBC00)
3446			return half(detail::binary, sign|0x3C00);
3447		switch(y.data_)
3448		{
3449			case 0x3800: return sqrt(x);
3450			case 0x3C00: return half(detail::binary, detail::check_underflow(x.data_));
3451			case 0x4000: return x * x;
3452			case 0xBC00: return half(detail::binary, 0x3C00) / x;
3453		}
3454		for(; absx<0x400; absx<<=1,--exp) ;
3455		detail::uint32 ilog = exp + (absx>>10), msign = detail::sign_mask(ilog), f, m = 
3456			(((ilog<<27)+((detail::log2(static_cast<detail::uint32>((absx&0x3FF)|0x400)<<20)+8)>>4))^msign) - msign;
3457		for(exp=-11; m<0x80000000; m<<=1,--exp) ;
3458		for(; absy<0x400; absy<<=1,--exp) ;
3459		m = detail::multiply64(m, static_cast<detail::uint32>((absy&0x3FF)|0x400)<<21);
3460		int i = m >> 31;
3461		exp += (absy>>10) + i;
3462		m <<= 1 - i;
3463		if(exp < 0)
3464		{
3465			f = m >> -exp;
3466			exp = 0;
3467		}
3468		else
3469		{
3470			f = (m<<exp) & 0x7FFFFFFF;
3471			exp = m >> (31-exp);
3472		}
3473		return half(detail::binary, detail::exp2_post<half::round_style>(f, exp, ((msign&1)^(y.data_>>15))!=0, sign));
3474	#endif
3475	}
3476
3477	/// \}
3478	/// \anchor trigonometric
3479	/// \name Trigonometric functions
3480	/// \{
3481
3482	/// Compute sine and cosine simultaneously.
3483	///	This returns the same results as sin() and cos() but is faster than calling each function individually.
3484	///
3485	/// This function is exact to rounding for all rounding modes.
3486	/// \param arg function argument
3487	/// \param sin variable to take sine of \a arg
3488	/// \param cos variable to take cosine of \a arg
3489	/// \exception FE_INVALID for signaling NaN or infinity
3490	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
3491	inline void sincos(half arg, half *sin, half *cos)
3492	{
3493	#ifdef HALF_ARITHMETIC_TYPE
3494		detail::internal_t f = detail::half2float<detail::internal_t>(arg.data_);
3495		*sin = half(detail::binary, detail::float2half<half::round_style>(std::sin(f)));
3496		*cos = half(detail::binary, detail::float2half<half::round_style>(std::cos(f)));
3497	#else
3498		int abs = arg.data_ & 0x7FFF, sign = arg.data_ >> 15, k;
3499		if(abs >= 0x7C00)
3500			*sin = *cos = half(detail::binary, (abs==0x7C00) ? detail::invalid() : detail::signal(arg.data_));
3501		else if(!abs)
3502		{
3503			*sin = arg;
3504			*cos = half(detail::binary, 0x3C00);
3505		}
3506		else if(abs < 0x2500)
3507		{
3508			*sin = half(detail::binary, detail::rounded<half::round_style,true>(arg.data_-1, 1, 1));
3509			*cos = half(detail::binary, detail::rounded<half::round_style,true>(0x3BFF, 1, 1));
3510		}
3511		else
3512		{
3513			if(half::round_style != std::round_to_nearest)
3514			{
3515				switch(abs)
3516				{
3517				case 0x48B7:
3518					*sin = half(detail::binary, detail::rounded<half::round_style,true>((~arg.data_&0x8000)|0x1D07, 1, 1));
3519					*cos = half(detail::binary, detail::rounded<half::round_style,true>(0xBBFF, 1, 1));
3520					return;
3521				case 0x598C:
3522					*sin = half(detail::binary, detail::rounded<half::round_style,true>((arg.data_&0x8000)|0x3BFF, 1, 1));
3523					*cos = half(detail::binary, detail::rounded<half::round_style,true>(0x80FC, 1, 1));
3524					return;
3525				case 0x6A64:
3526					*sin = half(detail::binary, detail::rounded<half::round_style,true>((~arg.data_&0x8000)|0x3BFE, 1, 1));
3527					*cos = half(detail::binary, detail::rounded<half::round_style,true>(0x27FF, 1, 1));
3528					return;
3529				case 0x6D8C:
3530					*sin = half(detail::binary, detail::rounded<half::round_style,true>((arg.data_&0x8000)|0x0FE6, 1, 1));
3531					*cos = half(detail::binary, detail::rounded<half::round_style,true>(0x3BFF, 1, 1));
3532					return;
3533				}
3534			}
3535			std::pair<detail::uint32,detail::uint32> sc = detail::sincos(detail::angle_arg(abs, k), 28);
3536			switch(k & 3)
3537			{
3538				case 1: sc = std::make_pair(sc.second, -sc.first); break;
3539				case 2: sc = std::make_pair(-sc.first, -sc.second); break;
3540				case 3: sc = std::make_pair(-sc.second, sc.first); break;
3541			}
3542			*sin = half(detail::binary, detail::fixed2half<half::round_style,30,true,true,true>((sc.first^-static_cast<detail::uint32>(sign))+sign));
3543			*cos = half(detail::binary, detail::fixed2half<half::round_style,30,true,true,true>(sc.second));
3544		}
3545	#endif
3546	}
3547
3548	/// Sine function.
3549	/// This function is exact to rounding for all rounding modes.
3550	///
3551	/// **See also:** Documentation for [std::sin](https://en.cppreference.com/w/cpp/numeric/math/sin).
3552	/// \param arg function argument
3553	/// \return sine value of \a arg
3554	/// \exception FE_INVALID for signaling NaN or infinity
3555	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
3556	inline half sin(half arg)
3557	{
3558	#ifdef HALF_ARITHMETIC_TYPE
3559		return half(detail::binary, detail::float2half<half::round_style>(std::sin(detail::half2float<detail::internal_t>(arg.data_))));
3560	#else
3561		int abs = arg.data_ & 0x7FFF, k;
3562		if(!abs)
3563			return arg;
3564		if(abs >= 0x7C00)
3565			return half(detail::binary, (abs==0x7C00) ? detail::invalid() : detail::signal(arg.data_));
3566		if(abs < 0x2900)
3567			return half(detail::binary, detail::rounded<half::round_style,true>(arg.data_-1, 1, 1));
3568		if(half::round_style != std::round_to_nearest)
3569			switch(abs)
3570			{
3571				case 0x48B7: return half(detail::binary, detail::rounded<half::round_style,true>((~arg.data_&0x8000)|0x1D07, 1, 1));
3572				case 0x6A64: return half(detail::binary, detail::rounded<half::round_style,true>((~arg.data_&0x8000)|0x3BFE, 1, 1));
3573				case 0x6D8C: return half(detail::binary, detail::rounded<half::round_style,true>((arg.data_&0x8000)|0x0FE6, 1, 1));
3574			}
3575		std::pair<detail::uint32,detail::uint32> sc = detail::sincos(detail::angle_arg(abs, k), 28);
3576		detail::uint32 sign = -static_cast<detail::uint32>(((k>>1)&1)^(arg.data_>>15));
3577		return half(detail::binary, detail::fixed2half<half::round_style,30,true,true,true>((((k&1) ? sc.second : sc.first)^sign) - sign));
3578	#endif
3579	}
3580
3581	/// Cosine function.
3582	/// This function is exact to rounding for all rounding modes.
3583	///
3584	/// **See also:** Documentation for [std::cos](https://en.cppreference.com/w/cpp/numeric/math/cos).
3585	/// \param arg function argument
3586	/// \return cosine value of \a arg
3587	/// \exception FE_INVALID for signaling NaN or infinity
3588	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
3589	inline half cos(half arg)
3590	{
3591	#ifdef HALF_ARITHMETIC_TYPE
3592		return half(detail::binary, detail::float2half<half::round_style>(std::cos(detail::half2float<detail::internal_t>(arg.data_))));
3593	#else
3594		int abs = arg.data_ & 0x7FFF, k;
3595		if(!abs)
3596			return half(detail::binary, 0x3C00);
3597		if(abs >= 0x7C00)
3598			return half(detail::binary, (abs==0x7C00) ? detail::invalid() : detail::signal(arg.data_));
3599		if(abs < 0x2500)
3600			return half(detail::binary, detail::rounded<half::round_style,true>(0x3BFF, 1, 1));
3601		if(half::round_style != std::round_to_nearest && abs == 0x598C)
3602			return half(detail::binary, detail::rounded<half::round_style,true>(0x80FC, 1, 1));
3603		std::pair<detail::uint32,detail::uint32> sc = detail::sincos(detail::angle_arg(abs, k), 28);
3604		detail::uint32 sign = -static_cast<detail::uint32>(((k>>1)^k)&1);
3605		return half(detail::binary, detail::fixed2half<half::round_style,30,true,true,true>((((k&1) ? sc.first : sc.second)^sign) - sign));
3606	#endif
3607	}
3608
3609	/// Tangent function.
3610	/// This function is exact to rounding for all rounding modes.
3611	///
3612	/// **See also:** Documentation for [std::tan](https://en.cppreference.com/w/cpp/numeric/math/tan).
3613	/// \param arg function argument
3614	/// \return tangent value of \a arg
3615	/// \exception FE_INVALID for signaling NaN or infinity
3616	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
3617	inline half tan(half arg)
3618	{
3619	#ifdef HALF_ARITHMETIC_TYPE
3620		return half(detail::binary, detail::float2half<half::round_style>(std::tan(detail::half2float<detail::internal_t>(arg.data_))));
3621	#else
3622		int abs = arg.data_ & 0x7FFF, exp = 13, k;
3623		if(!abs)
3624			return arg;
3625		if(abs >= 0x7C00)
3626			return half(detail::binary, (abs==0x7C00) ? detail::invalid() : detail::signal(arg.data_));
3627		if(abs < 0x2700)
3628			return half(detail::binary, detail::rounded<half::round_style,true>(arg.data_, 0, 1));
3629		if(half::round_style != std::round_to_nearest)
3630			switch(abs)
3631			{
3632				case 0x658C: return half(detail::binary, detail::rounded<half::round_style,true>((arg.data_&0x8000)|0x07E6, 1, 1));
3633				case 0x7330: return half(detail::binary, detail::rounded<half::round_style,true>((~arg.data_&0x8000)|0x4B62, 1, 1));
3634			}
3635		std::pair<detail::uint32,detail::uint32> sc = detail::sincos(detail::angle_arg(abs, k), 30);
3636		if(k & 1)
3637			sc = std::make_pair(-sc.second, sc.first);
3638		detail::uint32 signy = detail::sign_mask(sc.first), signx = detail::sign_mask(sc.second);
3639		detail::uint32 my = (sc.first^signy) - signy, mx = (sc.second^signx) - signx;
3640		for(; my<0x80000000; my<<=1,--exp) ;
3641		for(; mx<0x80000000; mx<<=1,++exp) ;
3642		return half(detail::binary, detail::tangent_post<half::round_style>(my, mx, exp, (signy^signx^arg.data_)&0x8000));
3643	#endif
3644	}
3645
3646	/// Arc sine.
3647	/// This function is exact to rounding for all rounding modes.
3648	///
3649	/// **See also:** Documentation for [std::asin](https://en.cppreference.com/w/cpp/numeric/math/asin).
3650	/// \param arg function argument
3651	/// \return arc sine value of \a arg
3652	/// \exception FE_INVALID for signaling NaN or if abs(\a arg) > 1
3653	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
3654	inline half asin(half arg)
3655	{
3656	#ifdef HALF_ARITHMETIC_TYPE
3657		return half(detail::binary, detail::float2half<half::round_style>(std::asin(detail::half2float<detail::internal_t>(arg.data_))));
3658	#else
3659		unsigned int abs = arg.data_ & 0x7FFF, sign = arg.data_ & 0x8000;
3660		if(!abs)
3661			return arg;
3662		if(abs >= 0x3C00)
3663			return half(detail::binary, (abs>0x7C00) ? detail::signal(arg.data_) : (abs>0x3C00) ? detail::invalid() :
3664										detail::rounded<half::round_style,true>(sign|0x3E48, 0, 1));
3665		if(abs < 0x2900)
3666			return half(detail::binary, detail::rounded<half::round_style,true>(arg.data_, 0, 1));
3667		if(half::round_style != std::round_to_nearest && (abs == 0x2B44 || abs == 0x2DC3))
3668			return half(detail::binary, detail::rounded<half::round_style,true>(arg.data_+1, 1, 1));
3669		std::pair<detail::uint32,detail::uint32> sc = detail::atan2_args(abs);
3670		detail::uint32 m = detail::atan2(sc.first, sc.second, (half::round_style==std::round_to_nearest) ? 27 : 26);
3671		return half(detail::binary, detail::fixed2half<half::round_style,30,false,true,true>(m, 14, sign));
3672	#endif
3673	}
3674
3675	/// Arc cosine function.
3676	/// This function is exact to rounding for all rounding modes.
3677	///
3678	/// **See also:** Documentation for [std::acos](https://en.cppreference.com/w/cpp/numeric/math/acos).
3679	/// \param arg function argument
3680	/// \return arc cosine value of \a arg
3681	/// \exception FE_INVALID for signaling NaN or if abs(\a arg) > 1
3682	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
3683	inline half acos(half arg)
3684	{
3685	#ifdef HALF_ARITHMETIC_TYPE
3686		return half(detail::binary, detail::float2half<half::round_style>(std::acos(detail::half2float<detail::internal_t>(arg.data_))));
3687	#else
3688		unsigned int abs = arg.data_ & 0x7FFF, sign = arg.data_ >> 15;
3689		if(!abs)
3690			return half(detail::binary, detail::rounded<half::round_style,true>(0x3E48, 0, 1));
3691		if(abs >= 0x3C00)
3692			return half(detail::binary,	(abs>0x7C00) ? detail::signal(arg.data_) : (abs>0x3C00) ? detail::invalid() :
3693										sign ? detail::rounded<half::round_style,true>(0x4248, 0, 1) : 0);
3694		std::pair<detail::uint32,detail::uint32> cs = detail::atan2_args(abs);
3695		detail::uint32 m = detail::atan2(cs.second, cs.first, 28);
3696		return half(detail::binary, detail::fixed2half<half::round_style,31,false,true,true>(sign ? (0xC90FDAA2-m) : m, 15, 0, sign));
3697	#endif
3698	}
3699
3700	/// Arc tangent function.
3701	/// This function is exact to rounding for all rounding modes.
3702	///
3703	/// **See also:** Documentation for [std::atan](https://en.cppreference.com/w/cpp/numeric/math/atan).
3704	/// \param arg function argument
3705	/// \return arc tangent value of \a arg
3706	/// \exception FE_INVALID for signaling NaN
3707	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
3708	inline half atan(half arg)
3709	{
3710	#ifdef HALF_ARITHMETIC_TYPE
3711		return half(detail::binary, detail::float2half<half::round_style>(std::atan(detail::half2float<detail::internal_t>(arg.data_))));
3712	#else
3713		unsigned int abs = arg.data_ & 0x7FFF, sign = arg.data_ & 0x8000;
3714		if(!abs)
3715			return arg;
3716		if(abs >= 0x7C00)
3717			return half(detail::binary, (abs==0x7C00) ? detail::rounded<half::round_style,true>(sign|0x3E48, 0, 1) : detail::signal(arg.data_));
3718		if(abs <= 0x2700)
3719			return half(detail::binary, detail::rounded<half::round_style,true>(arg.data_-1, 1, 1));
3720		int exp = (abs>>10) + (abs<=0x3FF);
3721		detail::uint32 my = (abs&0x3FF) | ((abs>0x3FF)<<10);
3722		detail::uint32 m = (exp>15) ?	detail::atan2(my<<19, 0x20000000>>(exp-15), (half::round_style==std::round_to_nearest) ? 26 : 24) :
3723										detail::atan2(my<<(exp+4), 0x20000000, (half::round_style==std::round_to_nearest) ? 30 : 28);
3724		return half(detail::binary, detail::fixed2half<half::round_style,30,false,true,true>(m, 14, sign));
3725	#endif
3726	}
3727
3728	/// Arc tangent function.
3729	/// This function may be 1 ULP off the correctly rounded exact result in ~0.005% of inputs for `std::round_to_nearest`, 
3730	/// in ~0.1% of inputs for `std::round_toward_zero` and in ~0.02% of inputs for any other rounding mode.
3731	///
3732	/// **See also:** Documentation for [std::atan2](https://en.cppreference.com/w/cpp/numeric/math/atan2).
3733	/// \param y numerator
3734	/// \param x denominator
3735	/// \return arc tangent value
3736	/// \exception FE_INVALID if \a x or \a y is signaling NaN
3737	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
3738	inline half atan2(half y, half x)
3739	{
3740	#ifdef HALF_ARITHMETIC_TYPE
3741		return half(detail::binary, detail::float2half<half::round_style>(std::atan2(detail::half2float<detail::internal_t>(y.data_), detail::half2float<detail::internal_t>(x.data_))));
3742	#else
3743		unsigned int absx = x.data_ & 0x7FFF, absy = y.data_ & 0x7FFF, signx = x.data_ >> 15, signy = y.data_ & 0x8000;
3744		if(absx >= 0x7C00 || absy >= 0x7C00)
3745		{
3746			if(absx > 0x7C00 || absy > 0x7C00)
3747				return half(detail::binary, detail::signal(x.data_, y.data_));
3748			if(absy == 0x7C00)
3749				return half(detail::binary, (absx<0x7C00) ?	detail::rounded<half::round_style,true>(signy|0x3E48, 0, 1) :
3750													signx ?	detail::rounded<half::round_style,true>(signy|0x40B6, 0, 1) :
3751															detail::rounded<half::round_style,true>(signy|0x3A48, 0, 1));
3752			return (x.data_==0x7C00) ? half(detail::binary, signy) : half(detail::binary, detail::rounded<half::round_style,true>(signy|0x4248, 0, 1));
3753		}
3754		if(!absy)
3755			return signx ? half(detail::binary, detail::rounded<half::round_style,true>(signy|0x4248, 0, 1)) : y;
3756		if(!absx)
3757			return half(detail::binary, detail::rounded<half::round_style,true>(signy|0x3E48, 0, 1));
3758		int d = (absy>>10) + (absy<=0x3FF) - (absx>>10) - (absx<=0x3FF);
3759		if(d > (signx ? 18 : 12))
3760			return half(detail::binary, detail::rounded<half::round_style,true>(signy|0x3E48, 0, 1));
3761		if(signx && d < -11)
3762			return half(detail::binary, detail::rounded<half::round_style,true>(signy|0x4248, 0, 1));
3763		if(!signx && d < ((half::round_style==std::round_toward_zero) ? -15 : -9))
3764		{
3765			for(; absy<0x400; absy<<=1,--d) ;
3766			detail::uint32 mx = ((absx<<1)&0x7FF) | 0x800, my = ((absy<<1)&0x7FF) | 0x800;
3767			int i = my < mx;
3768			d -= i;
3769			if(d < -25)
3770				return half(detail::binary, detail::underflow<half::round_style>(signy));
3771			my <<= 11 + i;
3772			return half(detail::binary, detail::fixed2half<half::round_style,11,false,false,true>(my/mx, d+14, signy, my%mx!=0));
3773		}
3774		detail::uint32 m = detail::atan2(	((absy&0x3FF)|((absy>0x3FF)<<10))<<(19+((d<0) ? d : (d>0) ? 0 : -1)),
3775											((absx&0x3FF)|((absx>0x3FF)<<10))<<(19-((d>0) ? d : (d<0) ? 0 : 1)));
3776		return half(detail::binary, detail::fixed2half<half::round_style,31,false,true,true>(signx ? (0xC90FDAA2-m) : m, 15, signy, signx));
3777	#endif
3778	}
3779
3780	/// \}
3781	/// \anchor hyperbolic
3782	/// \name Hyperbolic functions
3783	/// \{
3784
3785	/// Hyperbolic sine.
3786	/// This function is exact to rounding for all rounding modes.
3787	///
3788	/// **See also:** Documentation for [std::sinh](https://en.cppreference.com/w/cpp/numeric/math/sinh).
3789	/// \param arg function argument
3790	/// \return hyperbolic sine value of \a arg
3791	/// \exception FE_INVALID for signaling NaN
3792	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
3793	inline half sinh(half arg)
3794	{
3795	#ifdef HALF_ARITHMETIC_TYPE
3796		return half(detail::binary, detail::float2half<half::round_style>(std::sinh(detail::half2float<detail::internal_t>(arg.data_))));
3797	#else
3798		int abs = arg.data_ & 0x7FFF, exp;
3799		if(!abs || abs >= 0x7C00)
3800			return (abs>0x7C00) ? half(detail::binary, detail::signal(arg.data_)) : arg;
3801		if(abs <= 0x2900)
3802			return half(detail::binary, detail::rounded<half::round_style,true>(arg.data_, 0, 1));
3803		std::pair<detail::uint32,detail::uint32> mm = detail::hyperbolic_args(abs, exp, (half::round_style==std::round_to_nearest) ? 29 : 27);
3804		detail::uint32 m = mm.first - mm.second;
3805		for(exp+=13; m<0x80000000 && exp; m<<=1,--exp) ;
3806		unsigned int sign = arg.data_ & 0x8000;
3807		if(exp > 29)
3808			return half(detail::binary, detail::overflow<half::round_style>(sign));
3809		return half(detail::binary, detail::fixed2half<half::round_style,31,false,false,true>(m, exp, sign));
3810	#endif
3811	}
3812
3813	/// Hyperbolic cosine.
3814	/// This function is exact to rounding for all rounding modes.
3815	///
3816	/// **See also:** Documentation for [std::cosh](https://en.cppreference.com/w/cpp/numeric/math/cosh).
3817	/// \param arg function argument
3818	/// \return hyperbolic cosine value of \a arg
3819	/// \exception FE_INVALID for signaling NaN
3820	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
3821	inline half cosh(half arg)
3822	{
3823	#ifdef HALF_ARITHMETIC_TYPE
3824		return half(detail::binary, detail::float2half<half::round_style>(std::cosh(detail::half2float<detail::internal_t>(arg.data_))));
3825	#else
3826		int abs = arg.data_ & 0x7FFF, exp;
3827		if(!abs)
3828			return half(detail::binary, 0x3C00);
3829		if(abs >= 0x7C00)
3830			return half(detail::binary, (abs>0x7C00) ? detail::signal(arg.data_) : 0x7C00);
3831		std::pair<detail::uint32,detail::uint32> mm = detail::hyperbolic_args(abs, exp, (half::round_style==std::round_to_nearest) ? 23 : 26);
3832		detail::uint32 m = mm.first + mm.second, i = (~m&0xFFFFFFFF) >> 31;
3833		m = (m>>i) | (m&i) | 0x80000000;
3834		if((exp+=13+i) > 29)
3835			return half(detail::binary, detail::overflow<half::round_style>());
3836		return half(detail::binary, detail::fixed2half<half::round_style,31,false,false,true>(m, exp));
3837	#endif
3838	}
3839
3840	/// Hyperbolic tangent.
3841	/// This function is exact to rounding for all rounding modes.
3842	///
3843	/// **See also:** Documentation for [std::tanh](https://en.cppreference.com/w/cpp/numeric/math/tanh).
3844	/// \param arg function argument
3845	/// \return hyperbolic tangent value of \a arg
3846	/// \exception FE_INVALID for signaling NaN
3847	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
3848	inline half tanh(half arg)
3849	{
3850	#ifdef HALF_ARITHMETIC_TYPE
3851		return half(detail::binary, detail::float2half<half::round_style>(std::tanh(detail::half2float<detail::internal_t>(arg.data_))));
3852	#else
3853		int abs = arg.data_ & 0x7FFF, exp;
3854		if(!abs)
3855			return arg;
3856		if(abs >= 0x7C00)
3857			return half(detail::binary, (abs>0x7C00) ? detail::signal(arg.data_) : (arg.data_-0x4000));
3858		if(abs >= 0x4500)
3859			return half(detail::binary, detail::rounded<half::round_style,true>((arg.data_&0x8000)|0x3BFF, 1, 1));
3860		if(abs < 0x2700)
3861			return half(detail::binary, detail::rounded<half::round_style,true>(arg.data_-1, 1, 1));
3862		if(half::round_style != std::round_to_nearest && abs == 0x2D3F)
3863			return half(detail::binary, detail::rounded<half::round_style,true>(arg.data_-3, 0, 1));
3864		std::pair<detail::uint32,detail::uint32> mm = detail::hyperbolic_args(abs, exp, 27);
3865		detail::uint32 my = mm.first - mm.second - (half::round_style!=std::round_to_nearest), mx = mm.first + mm.second, i = (~mx&0xFFFFFFFF) >> 31;
3866		for(exp=13; my<0x80000000; my<<=1,--exp) ;
3867		mx = (mx>>i) | 0x80000000;
3868		return half(detail::binary, detail::tangent_post<half::round_style>(my, mx, exp-i, arg.data_&0x8000));
3869	#endif
3870	}
3871
3872	/// Hyperbolic area sine.
3873	/// This function is exact to rounding for all rounding modes.
3874	///
3875	/// **See also:** Documentation for [std::asinh](https://en.cppreference.com/w/cpp/numeric/math/asinh).
3876	/// \param arg function argument
3877	/// \return area sine value of \a arg
3878	/// \exception FE_INVALID for signaling NaN
3879	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
3880	inline half asinh(half arg)
3881	{
3882	#if defined(HALF_ARITHMETIC_TYPE) && HALF_ENABLE_CPP11_CMATH
3883		return half(detail::binary, detail::float2half<half::round_style>(std::asinh(detail::half2float<detail::internal_t>(arg.data_))));
3884	#else
3885		int abs = arg.data_ & 0x7FFF;
3886		if(!abs || abs >= 0x7C00)
3887			return (abs>0x7C00) ? half(detail::binary, detail::signal(arg.data_)) : arg;
3888		if(abs <= 0x2900)
3889			return half(detail::binary, detail::rounded<half::round_style,true>(arg.data_-1, 1, 1));
3890		if(half::round_style != std::round_to_nearest)
3891			switch(abs)
3892			{
3893				case 0x32D4: return half(detail::binary, detail::rounded<half::round_style,true>(arg.data_-13, 1, 1));
3894				case 0x3B5B: return half(detail::binary, detail::rounded<half::round_style,true>(arg.data_-197, 1, 1));
3895			}
3896		return half(detail::binary, detail::area<half::round_style,true>(arg.data_));
3897	#endif
3898	}
3899
3900	/// Hyperbolic area cosine.
3901	/// This function is exact to rounding for all rounding modes.
3902	///
3903	/// **See also:** Documentation for [std::acosh](https://en.cppreference.com/w/cpp/numeric/math/acosh).
3904	/// \param arg function argument
3905	/// \return area cosine value of \a arg
3906	/// \exception FE_INVALID for signaling NaN or arguments <1
3907	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
3908	inline half acosh(half arg)
3909	{
3910	#if defined(HALF_ARITHMETIC_TYPE) && HALF_ENABLE_CPP11_CMATH
3911		return half(detail::binary, detail::float2half<half::round_style>(std::acosh(detail::half2float<detail::internal_t>(arg.data_))));
3912	#else
3913		int abs = arg.data_ & 0x7FFF;
3914		if((arg.data_&0x8000) || abs < 0x3C00)
3915			return half(detail::binary, (abs<=0x7C00) ? detail::invalid() : detail::signal(arg.data_));
3916		if(abs == 0x3C00)
3917			return half(detail::binary, 0);
3918		if(arg.data_ >= 0x7C00)
3919			return (abs>0x7C00) ? half(detail::binary, detail::signal(arg.data_)) : arg;
3920		return half(detail::binary, detail::area<half::round_style,false>(arg.data_));
3921	#endif
3922	}
3923
3924	/// Hyperbolic area tangent.
3925	/// This function is exact to rounding for all rounding modes.
3926	///
3927	/// **See also:** Documentation for [std::atanh](https://en.cppreference.com/w/cpp/numeric/math/atanh).
3928	/// \param arg function argument
3929	/// \return area tangent value of \a arg
3930	/// \exception FE_INVALID for signaling NaN or if abs(\a arg) > 1
3931	/// \exception FE_DIVBYZERO for +/-1
3932	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
3933	inline half atanh(half arg)
3934	{
3935	#if defined(HALF_ARITHMETIC_TYPE) && HALF_ENABLE_CPP11_CMATH
3936		return half(detail::binary, detail::float2half<half::round_style>(std::atanh(detail::half2float<detail::internal_t>(arg.data_))));
3937	#else
3938		int abs = arg.data_ & 0x7FFF, exp = 0;
3939		if(!abs)
3940			return arg;
3941		if(abs >= 0x3C00)
3942			return half(detail::binary, (abs==0x3C00) ? detail::pole(arg.data_&0x8000) : (abs<=0x7C00) ? detail::invalid() : detail::signal(arg.data_));
3943		if(abs < 0x2700)
3944			return half(detail::binary, detail::rounded<half::round_style,true>(arg.data_, 0, 1));
3945		detail::uint32 m = static_cast<detail::uint32>((abs&0x3FF)|((abs>0x3FF)<<10)) << ((abs>>10)+(abs<=0x3FF)+6), my = 0x80000000 + m, mx = 0x80000000 - m;
3946		for(; mx<0x80000000; mx<<=1,++exp) ;
3947		int i = my >= mx, s;
3948		return half(detail::binary, detail::log2_post<half::round_style,0xB8AA3B2A>(detail::log2(
3949			(detail::divide64(my>>i, mx, s)+1)>>1, 27)+0x10, exp+i-1, 16, arg.data_&0x8000));
3950	#endif
3951	}
3952
3953	/// \}
3954	/// \anchor special
3955	/// \name Error and gamma functions
3956	/// \{
3957
3958	/// Error function.
3959	/// This function may be 1 ULP off the correctly rounded exact result for any rounding mode in <0.5% of inputs.
3960	///
3961	/// **See also:** Documentation for [std::erf](https://en.cppreference.com/w/cpp/numeric/math/erf).
3962	/// \param arg function argument
3963	/// \return error function value of \a arg
3964	/// \exception FE_INVALID for signaling NaN
3965	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
3966	inline half erf(half arg)
3967	{
3968	#if defined(HALF_ARITHMETIC_TYPE) && HALF_ENABLE_CPP11_CMATH
3969		return half(detail::binary, detail::float2half<half::round_style>(std::erf(detail::half2float<detail::internal_t>(arg.data_))));
3970	#else
3971		unsigned int abs = arg.data_ & 0x7FFF;
3972		if(!abs || abs >= 0x7C00)
3973			return (abs>=0x7C00) ? half(detail::binary, (abs==0x7C00) ? (arg.data_-0x4000) : detail::signal(arg.data_)) : arg;
3974		if(abs >= 0x4200)
3975			return half(detail::binary, detail::rounded<half::round_style,true>((arg.data_&0x8000)|0x3BFF, 1, 1));
3976		return half(detail::binary, detail::erf<half::round_style,false>(arg.data_));
3977	#endif
3978	}
3979
3980	/// Complementary error function.
3981	/// This function may be 1 ULP off the correctly rounded exact result for any rounding mode in <0.5% of inputs.
3982	///
3983	/// **See also:** Documentation for [std::erfc](https://en.cppreference.com/w/cpp/numeric/math/erfc).
3984	/// \param arg function argument
3985	/// \return 1 minus error function value of \a arg
3986	/// \exception FE_INVALID for signaling NaN
3987	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
3988	inline half erfc(half arg)
3989	{
3990	#if defined(HALF_ARITHMETIC_TYPE) && HALF_ENABLE_CPP11_CMATH
3991		return half(detail::binary, detail::float2half<half::round_style>(std::erfc(detail::half2float<detail::internal_t>(arg.data_))));
3992	#else
3993		unsigned int abs = arg.data_ & 0x7FFF, sign = arg.data_ & 0x8000;
3994		if(abs >= 0x7C00)
3995			return (abs>=0x7C00) ? half(detail::binary, (abs==0x7C00) ? (sign>>1) : detail::signal(arg.data_)) : arg;
3996		if(!abs)
3997			return half(detail::binary, 0x3C00);
3998		if(abs >= 0x4400)
3999			return half(detail::binary, detail::rounded<half::round_style,true>((sign>>1)-(sign>>15), sign>>15, 1));
4000		return half(detail::binary, detail::erf<half::round_style,true>(arg.data_));
4001	#endif
4002	}
4003
4004	/// Natural logarithm of gamma function.
4005	/// This function may be 1 ULP off the correctly rounded exact result for any rounding mode in ~0.025% of inputs.
4006	///
4007	/// **See also:** Documentation for [std::lgamma](https://en.cppreference.com/w/cpp/numeric/math/lgamma).
4008	/// \param arg function argument
4009	/// \return natural logarith of gamma function for \a arg
4010	/// \exception FE_INVALID for signaling NaN
4011	/// \exception FE_DIVBYZERO for 0 or negative integer arguments
4012	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
4013	inline half lgamma(half arg)
4014	{
4015	#if defined(HALF_ARITHMETIC_TYPE) && HALF_ENABLE_CPP11_CMATH
4016		return half(detail::binary, detail::float2half<half::round_style>(std::lgamma(detail::half2float<detail::internal_t>(arg.data_))));
4017	#else
4018		int abs = arg.data_ & 0x7FFF;
4019		if(abs >= 0x7C00)
4020			return half(detail::binary, (abs==0x7C00) ? 0x7C00 : detail::signal(arg.data_));
4021		if(!abs || arg.data_ >= 0xE400 || (arg.data_ >= 0xBC00 && !(abs&((1<<(25-(abs>>10)))-1))))
4022			return half(detail::binary, detail::pole());
4023		if(arg.data_ == 0x3C00 || arg.data_ == 0x4000)
4024			return half(detail::binary, 0);
4025		return half(detail::binary, detail::gamma<half::round_style,true>(arg.data_));
4026	#endif
4027	}
4028
4029	/// Gamma function.
4030	/// This function may be 1 ULP off the correctly rounded exact result for any rounding mode in <0.25% of inputs.
4031	///
4032	/// **See also:** Documentation for [std::tgamma](https://en.cppreference.com/w/cpp/numeric/math/tgamma).
4033	/// \param arg function argument
4034	/// \return gamma function value of \a arg
4035	/// \exception FE_INVALID for signaling NaN, negative infinity or negative integer arguments
4036	/// \exception FE_DIVBYZERO for 0
4037	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
4038	inline half tgamma(half arg)
4039	{
4040	#if defined(HALF_ARITHMETIC_TYPE) && HALF_ENABLE_CPP11_CMATH
4041		return half(detail::binary, detail::float2half<half::round_style>(std::tgamma(detail::half2float<detail::internal_t>(arg.data_))));
4042	#else
4043		unsigned int abs = arg.data_ & 0x7FFF;
4044		if(!abs)
4045			return half(detail::binary, detail::pole(arg.data_));
4046		if(abs >= 0x7C00)
4047			return (arg.data_==0x7C00) ? arg : half(detail::binary, detail::signal(arg.data_));
4048		if(arg.data_ >= 0xE400 || (arg.data_ >= 0xBC00 && !(abs&((1<<(25-(abs>>10)))-1))))
4049			return half(detail::binary, detail::invalid());
4050		if(arg.data_ >= 0xCA80)
4051			return half(detail::binary, detail::underflow<half::round_style>((1-((abs>>(25-(abs>>10)))&1))<<15));
4052		if(arg.data_ <= 0x100 || (arg.data_ >= 0x4900 && arg.data_ < 0x8000))
4053			return half(detail::binary, detail::overflow<half::round_style>());
4054		if(arg.data_ == 0x3C00)
4055			return arg;
4056		return half(detail::binary, detail::gamma<half::round_style,false>(arg.data_));
4057	#endif
4058	}
4059
4060	/// \}
4061	/// \anchor rounding
4062	/// \name Rounding
4063	/// \{
4064
4065	/// Nearest integer not less than half value.
4066	/// **See also:** Documentation for [std::ceil](https://en.cppreference.com/w/cpp/numeric/math/ceil).
4067	/// \param arg half to round
4068	/// \return nearest integer not less than \a arg
4069	/// \exception FE_INVALID for signaling NaN
4070	/// \exception FE_INEXACT if value had to be rounded
4071	inline half ceil(half arg) { return half(detail::binary, detail::integral<std::round_toward_infinity,true,true>(arg.data_)); }
4072
4073	/// Nearest integer not greater than half value.
4074	/// **See also:** Documentation for [std::floor](https://en.cppreference.com/w/cpp/numeric/math/floor).
4075	/// \param arg half to round
4076	/// \return nearest integer not greater than \a arg
4077	/// \exception FE_INVALID for signaling NaN
4078	/// \exception FE_INEXACT if value had to be rounded
4079	inline half floor(half arg) { return half(detail::binary, detail::integral<std::round_toward_neg_infinity,true,true>(arg.data_)); }
4080
4081	/// Nearest integer not greater in magnitude than half value.
4082	/// **See also:** Documentation for [std::trunc](https://en.cppreference.com/w/cpp/numeric/math/trunc).
4083	/// \param arg half to round
4084	/// \return nearest integer not greater in magnitude than \a arg
4085	/// \exception FE_INVALID for signaling NaN
4086	/// \exception FE_INEXACT if value had to be rounded
4087	inline half trunc(half arg) { return half(detail::binary, detail::integral<std::round_toward_zero,true,true>(arg.data_)); }
4088
4089	/// Nearest integer.
4090	/// **See also:** Documentation for [std::round](https://en.cppreference.com/w/cpp/numeric/math/round).
4091	/// \param arg half to round
4092	/// \return nearest integer, rounded away from zero in half-way cases
4093	/// \exception FE_INVALID for signaling NaN
4094	/// \exception FE_INEXACT if value had to be rounded
4095	inline half round(half arg) { return half(detail::binary, detail::integral<std::round_to_nearest,false,true>(arg.data_)); }
4096
4097	/// Nearest integer.
4098	/// **See also:** Documentation for [std::lround](https://en.cppreference.com/w/cpp/numeric/math/round).
4099	/// \param arg half to round
4100	/// \return nearest integer, rounded away from zero in half-way cases
4101	/// \exception FE_INVALID if value is not representable as `long`
4102	inline long lround(half arg) { return detail::half2int<std::round_to_nearest,false,false,long>(arg.data_); }
4103
4104	/// Nearest integer using half's internal rounding mode.
4105	/// **See also:** Documentation for [std::rint](https://en.cppreference.com/w/cpp/numeric/math/rint).
4106	/// \param arg half expression to round
4107	/// \return nearest integer using default rounding mode
4108	/// \exception FE_INVALID for signaling NaN
4109	/// \exception FE_INEXACT if value had to be rounded
4110	inline half rint(half arg) { return half(detail::binary, detail::integral<half::round_style,true,true>(arg.data_)); }
4111
4112	/// Nearest integer using half's internal rounding mode.
4113	/// **See also:** Documentation for [std::lrint](https://en.cppreference.com/w/cpp/numeric/math/rint).
4114	/// \param arg half expression to round
4115	/// \return nearest integer using default rounding mode
4116	/// \exception FE_INVALID if value is not representable as `long`
4117	/// \exception FE_INEXACT if value had to be rounded
4118	inline long lrint(half arg) { return detail::half2int<half::round_style,true,true,long>(arg.data_); }
4119
4120	/// Nearest integer using half's internal rounding mode.
4121	/// **See also:** Documentation for [std::nearbyint](https://en.cppreference.com/w/cpp/numeric/math/nearbyint).
4122	/// \param arg half expression to round
4123	/// \return nearest integer using default rounding mode
4124	/// \exception FE_INVALID for signaling NaN
4125	inline half nearbyint(half arg) { return half(detail::binary, detail::integral<half::round_style,true,false>(arg.data_)); }
4126#if HALF_ENABLE_CPP11_LONG_LONG
4127	/// Nearest integer.
4128	/// **See also:** Documentation for [std::llround](https://en.cppreference.com/w/cpp/numeric/math/round).
4129	/// \param arg half to round
4130	/// \return nearest integer, rounded away from zero in half-way cases
4131	/// \exception FE_INVALID if value is not representable as `long long`
4132	inline long long llround(half arg) { return detail::half2int<std::round_to_nearest,false,false,long long>(arg.data_); }
4133
4134	/// Nearest integer using half's internal rounding mode.
4135	/// **See also:** Documentation for [std::llrint](https://en.cppreference.com/w/cpp/numeric/math/rint).
4136	/// \param arg half expression to round
4137	/// \return nearest integer using default rounding mode
4138	/// \exception FE_INVALID if value is not representable as `long long`
4139	/// \exception FE_INEXACT if value had to be rounded
4140	inline long long llrint(half arg) { return detail::half2int<half::round_style,true,true,long long>(arg.data_); }
4141#endif
4142
4143	/// \}
4144	/// \anchor float
4145	/// \name Floating point manipulation
4146	/// \{
4147
4148	/// Decompress floating-point number.
4149	/// **See also:** Documentation for [std::frexp](https://en.cppreference.com/w/cpp/numeric/math/frexp).
4150	/// \param arg number to decompress
4151	/// \param exp address to store exponent at
4152	/// \return significant in range [0.5, 1)
4153	/// \exception FE_INVALID for signaling NaN
4154	inline half frexp(half arg, int *exp)
4155	{
4156		*exp = 0;
4157		unsigned int abs = arg.data_ & 0x7FFF;
4158		if(abs >= 0x7C00 || !abs)
4159			return (abs>0x7C00) ? half(detail::binary, detail::signal(arg.data_)) : arg;
4160		for(; abs<0x400; abs<<=1,--*exp) ;
4161		*exp += (abs>>10) - 14;
4162		return half(detail::binary, (arg.data_&0x8000)|0x3800|(abs&0x3FF));
4163	}
4164
4165	/// Multiply by power of two.
4166	/// This function is exact to rounding for all rounding modes.
4167	///
4168	/// **See also:** Documentation for [std::scalbln](https://en.cppreference.com/w/cpp/numeric/math/scalbn).
4169	/// \param arg number to modify
4170	/// \param exp power of two to multiply with
4171	/// \return \a arg multplied by 2 raised to \a exp
4172	/// \exception FE_INVALID for signaling NaN
4173	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
4174	inline half scalbln(half arg, long exp)
4175	{
4176		unsigned int abs = arg.data_ & 0x7FFF, sign = arg.data_ & 0x8000;
4177		if(abs >= 0x7C00 || !abs)
4178			return (abs>0x7C00) ? half(detail::binary, detail::signal(arg.data_)) : arg;
4179		for(; abs<0x400; abs<<=1,--exp) ;
4180		exp += abs >> 10;
4181		if(exp > 30)
4182			return half(detail::binary, detail::overflow<half::round_style>(sign));
4183		else if(exp < -10)
4184			return half(detail::binary, detail::underflow<half::round_style>(sign));
4185		else if(exp > 0)
4186			return half(detail::binary, sign|(exp<<10)|(abs&0x3FF));
4187		unsigned int m = (abs&0x3FF) | 0x400;
4188		return half(detail::binary, detail::rounded<half::round_style,false>(sign|(m>>(1-exp)), (m>>-exp)&1, (m&((1<<-exp)-1))!=0));
4189	}
4190
4191	/// Multiply by power of two.
4192	/// This function is exact to rounding for all rounding modes.
4193	///
4194	/// **See also:** Documentation for [std::scalbn](https://en.cppreference.com/w/cpp/numeric/math/scalbn).
4195	/// \param arg number to modify
4196	/// \param exp power of two to multiply with
4197	/// \return \a arg multplied by 2 raised to \a exp
4198	/// \exception FE_INVALID for signaling NaN
4199	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
4200	inline half scalbn(half arg, int exp) { return scalbln(arg, exp); }
4201
4202	/// Multiply by power of two.
4203	/// This function is exact to rounding for all rounding modes.
4204	///
4205	/// **See also:** Documentation for [std::ldexp](https://en.cppreference.com/w/cpp/numeric/math/ldexp).
4206	/// \param arg number to modify
4207	/// \param exp power of two to multiply with
4208	/// \return \a arg multplied by 2 raised to \a exp
4209	/// \exception FE_INVALID for signaling NaN
4210	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
4211	inline half ldexp(half arg, int exp) { return scalbln(arg, exp); }
4212
4213	/// Extract integer and fractional parts.
4214	/// **See also:** Documentation for [std::modf](https://en.cppreference.com/w/cpp/numeric/math/modf).
4215	/// \param arg number to decompress
4216	/// \param iptr address to store integer part at
4217	/// \return fractional part
4218	/// \exception FE_INVALID for signaling NaN
4219	inline half modf(half arg, half *iptr)
4220	{
4221		unsigned int abs = arg.data_ & 0x7FFF;
4222		if(abs > 0x7C00)
4223		{
4224			arg = half(detail::binary, detail::signal(arg.data_));
4225			return *iptr = arg, arg;
4226		}
4227		if(abs >= 0x6400)
4228			return *iptr = arg, half(detail::binary, arg.data_&0x8000);
4229		if(abs < 0x3C00)
4230			return iptr->data_ = arg.data_ & 0x8000, arg;
4231		unsigned int exp = abs >> 10, mask = (1<<(25-exp)) - 1, m = arg.data_ & mask;
4232		iptr->data_ = arg.data_ & ~mask;
4233		if(!m)
4234			return half(detail::binary, arg.data_&0x8000);
4235		for(; m<0x400; m<<=1,--exp) ;
4236		return half(detail::binary, (arg.data_&0x8000)|(exp<<10)|(m&0x3FF));
4237	}
4238
4239	/// Extract exponent.
4240	/// **See also:** Documentation for [std::ilogb](https://en.cppreference.com/w/cpp/numeric/math/ilogb).
4241	/// \param arg number to query
4242	/// \return floating-point exponent
4243	/// \retval FP_ILOGB0 for zero
4244	/// \retval FP_ILOGBNAN for NaN
4245	/// \retval INT_MAX for infinity
4246	/// \exception FE_INVALID for 0 or infinite values
4247	inline int ilogb(half arg)
4248	{
4249		int abs = arg.data_ & 0x7FFF, exp;
4250		if(!abs || abs >= 0x7C00)
4251		{
4252			detail::raise(FE_INVALID);
4253			return !abs ? FP_ILOGB0 : (abs==0x7C00) ? INT_MAX : FP_ILOGBNAN;
4254		}
4255		for(exp=(abs>>10)-15; abs<0x200; abs<<=1,--exp) ;
4256		return exp;
4257	}
4258
4259	/// Extract exponent.
4260	/// **See also:** Documentation for [std::logb](https://en.cppreference.com/w/cpp/numeric/math/logb).
4261	/// \param arg number to query
4262	/// \return floating-point exponent
4263	/// \exception FE_INVALID for signaling NaN
4264	/// \exception FE_DIVBYZERO for 0
4265	inline half logb(half arg)
4266	{
4267		int abs = arg.data_ & 0x7FFF, exp;
4268		if(!abs)
4269			return half(detail::binary, detail::pole(0x8000));
4270		if(abs >= 0x7C00)
4271			return half(detail::binary, (abs==0x7C00) ? 0x7C00 : detail::signal(arg.data_));
4272		for(exp=(abs>>10)-15; abs<0x200; abs<<=1,--exp) ;
4273		unsigned int value = static_cast<unsigned>(exp<0) << 15;
4274		if(exp)
4275		{
4276			unsigned int m = std::abs(exp) << 6;
4277			for(exp=18; m<0x400; m<<=1,--exp) ;
4278			value |= (exp<<10) + m;
4279		}
4280		return half(detail::binary, value);
4281	}
4282
4283	/// Next representable value.
4284	/// **See also:** Documentation for [std::nextafter](https://en.cppreference.com/w/cpp/numeric/math/nextafter).
4285	/// \param from value to compute next representable value for
4286	/// \param to direction towards which to compute next value
4287	/// \return next representable value after \a from in direction towards \a to
4288	/// \exception FE_INVALID for signaling NaN
4289	/// \exception FE_OVERFLOW for infinite result from finite argument
4290	/// \exception FE_UNDERFLOW for subnormal result
4291	inline half nextafter(half from, half to)
4292	{
4293		int fabs = from.data_ & 0x7FFF, tabs = to.data_ & 0x7FFF;
4294		if(fabs > 0x7C00 || tabs > 0x7C00)
4295			return half(detail::binary, detail::signal(from.data_, to.data_));
4296		if(from.data_ == to.data_ || !(fabs|tabs))
4297			return to;
4298		if(!fabs)
4299		{
4300			detail::raise(FE_UNDERFLOW, !HALF_ERRHANDLING_UNDERFLOW_TO_INEXACT);
4301			return half(detail::binary, (to.data_&0x8000)+1);
4302		}
4303		unsigned int out = from.data_ + (((from.data_>>15)^static_cast<unsigned>(
4304			(from.data_^(0x8000|(0x8000-(from.data_>>15))))<(to.data_^(0x8000|(0x8000-(to.data_>>15))))))<<1) - 1;
4305		detail::raise(FE_OVERFLOW, fabs<0x7C00 && (out&0x7C00)==0x7C00);
4306		detail::raise(FE_UNDERFLOW, !HALF_ERRHANDLING_UNDERFLOW_TO_INEXACT && (out&0x7C00)<0x400);
4307		return half(detail::binary, out);
4308	}
4309
4310	/// Next representable value.
4311	/// **See also:** Documentation for [std::nexttoward](https://en.cppreference.com/w/cpp/numeric/math/nexttoward).
4312	/// \param from value to compute next representable value for
4313	/// \param to direction towards which to compute next value
4314	/// \return next representable value after \a from in direction towards \a to
4315	/// \exception FE_INVALID for signaling NaN
4316	/// \exception FE_OVERFLOW for infinite result from finite argument
4317	/// \exception FE_UNDERFLOW for subnormal result
4318	inline half nexttoward(half from, long double to)
4319	{
4320		int fabs = from.data_ & 0x7FFF;
4321		if(fabs > 0x7C00)
4322			return half(detail::binary, detail::signal(from.data_));
4323		long double lfrom = static_cast<long double>(from);
4324		if(detail::builtin_isnan(to) || lfrom == to)
4325			return half(static_cast<float>(to));
4326		if(!fabs)
4327		{
4328			detail::raise(FE_UNDERFLOW, !HALF_ERRHANDLING_UNDERFLOW_TO_INEXACT);
4329			return half(detail::binary, (static_cast<unsigned>(detail::builtin_signbit(to))<<15)+1);
4330		}
4331		unsigned int out = from.data_ + (((from.data_>>15)^static_cast<unsigned>(lfrom<to))<<1) - 1;
4332		detail::raise(FE_OVERFLOW, (out&0x7FFF)==0x7C00);
4333		detail::raise(FE_UNDERFLOW, !HALF_ERRHANDLING_UNDERFLOW_TO_INEXACT && (out&0x7FFF)<0x400);
4334		return half(detail::binary, out);
4335	}
4336
4337	/// Take sign.
4338	/// **See also:** Documentation for [std::copysign](https://en.cppreference.com/w/cpp/numeric/math/copysign).
4339	/// \param x value to change sign for
4340	/// \param y value to take sign from
4341	/// \return value equal to \a x in magnitude and to \a y in sign
4342	inline HALF_CONSTEXPR half copysign(half x, half y) { return half(detail::binary, x.data_^((x.data_^y.data_)&0x8000)); }
4343
4344	/// \}
4345	/// \anchor classification
4346	/// \name Floating point classification
4347	/// \{
4348
4349	/// Classify floating-point value.
4350	/// **See also:** Documentation for [std::fpclassify](https://en.cppreference.com/w/cpp/numeric/math/fpclassify).
4351	/// \param arg number to classify
4352	/// \retval FP_ZERO for positive and negative zero
4353	/// \retval FP_SUBNORMAL for subnormal numbers
4354	/// \retval FP_INFINITY for positive and negative infinity
4355	/// \retval FP_NAN for NaNs
4356	/// \retval FP_NORMAL for all other (normal) values
4357	inline HALF_CONSTEXPR int fpclassify(half arg)
4358	{
4359		return	!(arg.data_&0x7FFF) ? FP_ZERO :
4360				((arg.data_&0x7FFF)<0x400) ? FP_SUBNORMAL :
4361				((arg.data_&0x7FFF)<0x7C00) ? FP_NORMAL :
4362				((arg.data_&0x7FFF)==0x7C00) ? FP_INFINITE :
4363				FP_NAN;
4364	}
4365
4366	/// Check if finite number.
4367	/// **See also:** Documentation for [std::isfinite](https://en.cppreference.com/w/cpp/numeric/math/isfinite).
4368	/// \param arg number to check
4369	/// \retval true if neither infinity nor NaN
4370	/// \retval false else
4371	inline HALF_CONSTEXPR bool isfinite(half arg) { return (arg.data_&0x7C00) != 0x7C00; }
4372
4373	/// Check for infinity.
4374	/// **See also:** Documentation for [std::isinf](https://en.cppreference.com/w/cpp/numeric/math/isinf).
4375	/// \param arg number to check
4376	/// \retval true for positive or negative infinity
4377	/// \retval false else
4378	inline HALF_CONSTEXPR bool isinf(half arg) { return (arg.data_&0x7FFF) == 0x7C00; }
4379
4380	/// Check for NaN.
4381	/// **See also:** Documentation for [std::isnan](https://en.cppreference.com/w/cpp/numeric/math/isnan).
4382	/// \param arg number to check
4383	/// \retval true for NaNs
4384	/// \retval false else
4385	inline HALF_CONSTEXPR bool isnan(half arg) { return (arg.data_&0x7FFF) > 0x7C00; }
4386
4387	/// Check if normal number.
4388	/// **See also:** Documentation for [std::isnormal](https://en.cppreference.com/w/cpp/numeric/math/isnormal).
4389	/// \param arg number to check
4390	/// \retval true if normal number
4391	/// \retval false if either subnormal, zero, infinity or NaN
4392	inline HALF_CONSTEXPR bool isnormal(half arg) { return ((arg.data_&0x7C00)!=0) & ((arg.data_&0x7C00)!=0x7C00); }
4393
4394	/// Check sign.
4395	/// **See also:** Documentation for [std::signbit](https://en.cppreference.com/w/cpp/numeric/math/signbit).
4396	/// \param arg number to check
4397	/// \retval true for negative number
4398	/// \retval false for positive number
4399	inline HALF_CONSTEXPR bool signbit(half arg) { return (arg.data_&0x8000) != 0; }
4400
4401	/// \}
4402	/// \anchor compfunc
4403	/// \name Comparison
4404	/// \{
4405
4406	/// Quiet comparison for greater than.
4407	/// **See also:** Documentation for [std::isgreater](https://en.cppreference.com/w/cpp/numeric/math/isgreater).
4408	/// \param x first operand
4409	/// \param y second operand
4410	/// \retval true if \a x greater than \a y
4411	/// \retval false else
4412	inline HALF_CONSTEXPR bool isgreater(half x, half y)
4413	{
4414		return ((x.data_^(0x8000|(0x8000-(x.data_>>15))))+(x.data_>>15)) > ((y.data_^(0x8000|(0x8000-(y.data_>>15))))+(y.data_>>15)) && !isnan(x) && !isnan(y);
4415	}
4416
4417	/// Quiet comparison for greater equal.
4418	/// **See also:** Documentation for [std::isgreaterequal](https://en.cppreference.com/w/cpp/numeric/math/isgreaterequal).
4419	/// \param x first operand
4420	/// \param y second operand
4421	/// \retval true if \a x greater equal \a y
4422	/// \retval false else
4423	inline HALF_CONSTEXPR bool isgreaterequal(half x, half y)
4424	{
4425		return ((x.data_^(0x8000|(0x8000-(x.data_>>15))))+(x.data_>>15)) >= ((y.data_^(0x8000|(0x8000-(y.data_>>15))))+(y.data_>>15)) && !isnan(x) && !isnan(y);
4426	}
4427
4428	/// Quiet comparison for less than.
4429	/// **See also:** Documentation for [std::isless](https://en.cppreference.com/w/cpp/numeric/math/isless).
4430	/// \param x first operand
4431	/// \param y second operand
4432	/// \retval true if \a x less than \a y
4433	/// \retval false else
4434	inline HALF_CONSTEXPR bool isless(half x, half y)
4435	{
4436		return ((x.data_^(0x8000|(0x8000-(x.data_>>15))))+(x.data_>>15)) < ((y.data_^(0x8000|(0x8000-(y.data_>>15))))+(y.data_>>15)) && !isnan(x) && !isnan(y);
4437	}
4438
4439	/// Quiet comparison for less equal.
4440	/// **See also:** Documentation for [std::islessequal](https://en.cppreference.com/w/cpp/numeric/math/islessequal).
4441	/// \param x first operand
4442	/// \param y second operand
4443	/// \retval true if \a x less equal \a y
4444	/// \retval false else
4445	inline HALF_CONSTEXPR bool islessequal(half x, half y)
4446	{
4447		return ((x.data_^(0x8000|(0x8000-(x.data_>>15))))+(x.data_>>15)) <= ((y.data_^(0x8000|(0x8000-(y.data_>>15))))+(y.data_>>15)) && !isnan(x) && !isnan(y);
4448	}
4449
4450	/// Quiet comarison for less or greater.
4451	/// **See also:** Documentation for [std::islessgreater](https://en.cppreference.com/w/cpp/numeric/math/islessgreater).
4452	/// \param x first operand
4453	/// \param y second operand
4454	/// \retval true if either less or greater
4455	/// \retval false else
4456	inline HALF_CONSTEXPR bool islessgreater(half x, half y)
4457	{
4458		return x.data_!=y.data_ && ((x.data_|y.data_)&0x7FFF) && !isnan(x) && !isnan(y);
4459	}
4460
4461	/// Quiet check if unordered.
4462	/// **See also:** Documentation for [std::isunordered](https://en.cppreference.com/w/cpp/numeric/math/isunordered).
4463	/// \param x first operand
4464	/// \param y second operand
4465	/// \retval true if unordered (one or two NaN operands)
4466	/// \retval false else
4467	inline HALF_CONSTEXPR bool isunordered(half x, half y) { return isnan(x) || isnan(y); }
4468
4469	/// \}
4470	/// \anchor casting
4471	/// \name Casting
4472	/// \{
4473
4474	/// Cast to or from half-precision floating-point number.
4475	/// This casts between [half](\ref half_float::half) and any built-in arithmetic type. The values are converted 
4476	/// directly using the default rounding mode, without any roundtrip over `float` that a `static_cast` would otherwise do.
4477	///
4478	/// Using this cast with neither of the two types being a [half](\ref half_float::half) or with any of the two types 
4479	/// not being a built-in arithmetic type (apart from [half](\ref half_float::half), of course) results in a compiler 
4480	/// error and casting between [half](\ref half_float::half)s returns the argument unmodified.
4481	/// \tparam T destination type (half or built-in arithmetic type)
4482	/// \tparam U source type (half or built-in arithmetic type)
4483	/// \param arg value to cast
4484	/// \return \a arg converted to destination type
4485	/// \exception FE_INVALID if \a T is integer type and result is not representable as \a T
4486	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
4487	template<typename T,typename U> T half_cast(U arg) { return detail::half_caster<T,U>::cast(arg); }
4488
4489	/// Cast to or from half-precision floating-point number.
4490	/// This casts between [half](\ref half_float::half) and any built-in arithmetic type. The values are converted 
4491	/// directly using the specified rounding mode, without any roundtrip over `float` that a `static_cast` would otherwise do.
4492	///
4493	/// Using this cast with neither of the two types being a [half](\ref half_float::half) or with any of the two types 
4494	/// not being a built-in arithmetic type (apart from [half](\ref half_float::half), of course) results in a compiler 
4495	/// error and casting between [half](\ref half_float::half)s returns the argument unmodified.
4496	/// \tparam T destination type (half or built-in arithmetic type)
4497	/// \tparam R rounding mode to use.
4498	/// \tparam U source type (half or built-in arithmetic type)
4499	/// \param arg value to cast
4500	/// \return \a arg converted to destination type
4501	/// \exception FE_INVALID if \a T is integer type and result is not representable as \a T
4502	/// \exception FE_OVERFLOW, ...UNDERFLOW, ...INEXACT according to rounding
4503	template<typename T,std::float_round_style R,typename U> T half_cast(U arg) { return detail::half_caster<T,U,R>::cast(arg); }
4504	/// \}
4505
4506	/// \}
4507	/// \anchor errors
4508	/// \name Error handling
4509	/// \{
4510
4511	/// Clear exception flags.
4512	/// This function works even if [automatic exception flag handling](\ref HALF_ERRHANDLING_FLAGS) is disabled, 
4513	/// but in that case manual flag management is the only way to raise flags.
4514	///
4515	/// **See also:** Documentation for [std::feclearexcept](https://en.cppreference.com/w/cpp/numeric/fenv/feclearexcept).
4516	/// \param excepts OR of exceptions to clear
4517	/// \retval 0 all selected flags cleared successfully
4518	inline int feclearexcept(int excepts) { detail::errflags() &= ~excepts; return 0; }
4519
4520	/// Test exception flags.
4521	/// This function works even if [automatic exception flag handling](\ref HALF_ERRHANDLING_FLAGS) is disabled, 
4522	/// but in that case manual flag management is the only way to raise flags.
4523	///
4524	/// **See also:** Documentation for [std::fetestexcept](https://en.cppreference.com/w/cpp/numeric/fenv/fetestexcept).
4525	/// \param excepts OR of exceptions to test
4526	/// \return OR of selected exceptions if raised
4527	inline int fetestexcept(int excepts) { return detail::errflags() & excepts; }
4528
4529	/// Raise exception flags.
4530	/// This raises the specified floating point exceptions and also invokes any additional automatic exception handling as 
4531	/// configured with the [HALF_ERRHANDLIG_...](\ref HALF_ERRHANDLING_ERRNO) preprocessor symbols.
4532	/// This function works even if [automatic exception flag handling](\ref HALF_ERRHANDLING_FLAGS) is disabled, 
4533	/// but in that case manual flag management is the only way to raise flags.
4534	///
4535	/// **See also:** Documentation for [std::feraiseexcept](https://en.cppreference.com/w/cpp/numeric/fenv/feraiseexcept).
4536	/// \param excepts OR of exceptions to raise
4537	/// \retval 0 all selected exceptions raised successfully
4538	inline int feraiseexcept(int excepts) { detail::errflags() |= excepts; detail::raise(excepts); return 0; }
4539
4540	/// Save exception flags.
4541	/// This function works even if [automatic exception flag handling](\ref HALF_ERRHANDLING_FLAGS) is disabled, 
4542	/// but in that case manual flag management is the only way to raise flags.
4543	///
4544	/// **See also:** Documentation for [std::fegetexceptflag](https://en.cppreference.com/w/cpp/numeric/fenv/feexceptflag).
4545	/// \param flagp adress to store flag state at
4546	/// \param excepts OR of flags to save
4547	/// \retval 0 for success
4548	inline int fegetexceptflag(int *flagp, int excepts) { *flagp = detail::errflags() & excepts; return 0; }
4549
4550	/// Restore exception flags.
4551	/// This only copies the specified exception state (including unset flags) without incurring any additional exception handling.
4552	/// This function works even if [automatic exception flag handling](\ref HALF_ERRHANDLING_FLAGS) is disabled, 
4553	/// but in that case manual flag management is the only way to raise flags.
4554	///
4555	/// **See also:** Documentation for [std::fesetexceptflag](https://en.cppreference.com/w/cpp/numeric/fenv/feexceptflag).
4556	/// \param flagp adress to take flag state from
4557	/// \param excepts OR of flags to restore
4558	/// \retval 0 for success
4559	inline int fesetexceptflag(const int *flagp, int excepts) { detail::errflags() = (detail::errflags()|(*flagp&excepts)) & (*flagp|~excepts); return 0; }
4560
4561	/// Throw C++ exceptions based on set exception flags.
4562	/// This function manually throws a corresponding C++ exception if one of the specified flags is set, 
4563	/// no matter if automatic throwing (via [HALF_ERRHANDLING_THROW_...](\ref HALF_ERRHANDLING_THROW_INVALID)) is enabled or not.
4564	/// This function works even if [automatic exception flag handling](\ref HALF_ERRHANDLING_FLAGS) is disabled, 
4565	/// but in that case manual flag management is the only way to raise flags.
4566	/// \param excepts OR of exceptions to test
4567	/// \param msg error message to use for exception description
4568	/// \throw std::domain_error if `FE_INVALID` or `FE_DIVBYZERO` is selected and set
4569	/// \throw std::overflow_error if `FE_OVERFLOW` is selected and set
4570	/// \throw std::underflow_error if `FE_UNDERFLOW` is selected and set
4571	/// \throw std::range_error if `FE_INEXACT` is selected and set
4572	inline void fethrowexcept(int excepts, const char *msg = "")
4573	{
4574		excepts &= detail::errflags();
4575		if(excepts & (FE_INVALID|FE_DIVBYZERO))
4576			throw std::domain_error(msg);
4577		if(excepts & FE_OVERFLOW)
4578			throw std::overflow_error(msg);
4579		if(excepts & FE_UNDERFLOW)
4580			throw std::underflow_error(msg);
4581		if(excepts & FE_INEXACT)
4582			throw std::range_error(msg);
4583	}
4584	/// \}
4585}
4586
4587
4588#undef HALF_UNUSED_NOERR
4589#undef HALF_CONSTEXPR
4590#undef HALF_CONSTEXPR_CONST
4591#undef HALF_CONSTEXPR_NOERR
4592#undef HALF_NOEXCEPT
4593#undef HALF_NOTHROW
4594#undef HALF_THREAD_LOCAL
4595#undef HALF_TWOS_COMPLEMENT_INT
4596#ifdef HALF_POP_WARNINGS
4597	#pragma warning(pop)
4598	#undef HALF_POP_WARNINGS
4599#endif
4600
4601#endif