cpp_dec_float.hpp 130 KB

12345678910111213141516171819202122232425262728293031323334353637383940414243444546474849505152535455565758596061626364656667686970717273747576777879808182838485868788899091929394959697989910010110210310410510610710810911011111211311411511611711811912012112212312412512612712812913013113213313413513613713813914014114214314414514614714814915015115215315415515615715815916016116216316416516616716816917017117217317417517617717817918018118218318418518618718818919019119219319419519619719819920020120220320420520620720820921021121221321421521621721821922022122222322422522622722822923023123223323423523623723823924024124224324424524624724824925025125225325425525625725825926026126226326426526626726826927027127227327427527627727827928028128228328428528628728828929029129229329429529629729829930030130230330430530630730830931031131231331431531631731831932032132232332432532632732832933033133233333433533633733833934034134234334434534634734834935035135235335435535635735835936036136236336436536636736836937037137237337437537637737837938038138238338438538638738838939039139239339439539639739839940040140240340440540640740840941041141241341441541641741841942042142242342442542642742842943043143243343443543643743843944044144244344444544644744844945045145245345445545645745845946046146246346446546646746846947047147247347447547647747847948048148248348448548648748848949049149249349449549649749849950050150250350450550650750850951051151251351451551651751851952052152252352452552652752852953053153253353453553653753853954054154254354454554654754854955055155255355455555655755855956056156256356456556656756856957057157257357457557657757857958058158258358458558658758858959059159259359459559659759859960060160260360460560660760860961061161261361461561661761861962062162262362462562662762862963063163263363463563663763863964064164264364464564664764864965065165265365465565665765865966066166266366466566666766866967067167267367467567667767867968068168268368468568668768868969069169269369469569669769869970070170270370470570670770870971071171271371471571671771871972072172272372472572672772872973073173273373473573673773873974074174274374474574674774874975075175275375475575675775875976076176276376476576676776876977077177277377477577677777877978078178278378478578678778878979079179279379479579679779879980080180280380480580680780880981081181281381481581681781881982082182282382482582682782882983083183283383483583683783883984084184284384484584684784884985085185285385485585685785885986086186286386486586686786886987087187287387487587687787887988088188288388488588688788888989089189289389489589689789889990090190290390490590690790890991091191291391491591691791891992092192292392492592692792892993093193293393493593693793893994094194294394494594694794894995095195295395495595695795895996096196296396496596696796896997097197297397497597697797897998098198298398498598698798898999099199299399499599699799899910001001100210031004100510061007100810091010101110121013101410151016101710181019102010211022102310241025102610271028102910301031103210331034103510361037103810391040104110421043104410451046104710481049105010511052105310541055105610571058105910601061106210631064106510661067106810691070107110721073107410751076107710781079108010811082108310841085108610871088108910901091109210931094109510961097109810991100110111021103110411051106110711081109111011111112111311141115111611171118111911201121112211231124112511261127112811291130113111321133113411351136113711381139114011411142114311441145114611471148114911501151115211531154115511561157115811591160116111621163116411651166116711681169117011711172117311741175117611771178117911801181118211831184118511861187118811891190119111921193119411951196119711981199120012011202120312041205120612071208120912101211121212131214121512161217121812191220122112221223122412251226122712281229123012311232123312341235123612371238123912401241124212431244124512461247124812491250125112521253125412551256125712581259126012611262126312641265126612671268126912701271127212731274127512761277127812791280128112821283128412851286128712881289129012911292129312941295129612971298129913001301130213031304130513061307130813091310131113121313131413151316131713181319132013211322132313241325132613271328132913301331133213331334133513361337133813391340134113421343134413451346134713481349135013511352135313541355135613571358135913601361136213631364136513661367136813691370137113721373137413751376137713781379138013811382138313841385138613871388138913901391139213931394139513961397139813991400140114021403140414051406140714081409141014111412141314141415141614171418141914201421142214231424142514261427142814291430143114321433143414351436143714381439144014411442144314441445144614471448144914501451145214531454145514561457145814591460146114621463146414651466146714681469147014711472147314741475147614771478147914801481148214831484148514861487148814891490149114921493149414951496149714981499150015011502150315041505150615071508150915101511151215131514151515161517151815191520152115221523152415251526152715281529153015311532153315341535153615371538153915401541154215431544154515461547154815491550155115521553155415551556155715581559156015611562156315641565156615671568156915701571157215731574157515761577157815791580158115821583158415851586158715881589159015911592159315941595159615971598159916001601160216031604160516061607160816091610161116121613161416151616161716181619162016211622162316241625162616271628162916301631163216331634163516361637163816391640164116421643164416451646164716481649165016511652165316541655165616571658165916601661166216631664166516661667166816691670167116721673167416751676167716781679168016811682168316841685168616871688168916901691169216931694169516961697169816991700170117021703170417051706170717081709171017111712171317141715171617171718171917201721172217231724172517261727172817291730173117321733173417351736173717381739174017411742174317441745174617471748174917501751175217531754175517561757175817591760176117621763176417651766176717681769177017711772177317741775177617771778177917801781178217831784178517861787178817891790179117921793179417951796179717981799180018011802180318041805180618071808180918101811181218131814181518161817181818191820182118221823182418251826182718281829183018311832183318341835183618371838183918401841184218431844184518461847184818491850185118521853185418551856185718581859186018611862186318641865186618671868186918701871187218731874187518761877187818791880188118821883188418851886188718881889189018911892189318941895189618971898189919001901190219031904190519061907190819091910191119121913191419151916191719181919192019211922192319241925192619271928192919301931193219331934193519361937193819391940194119421943194419451946194719481949195019511952195319541955195619571958195919601961196219631964196519661967196819691970197119721973197419751976197719781979198019811982198319841985198619871988198919901991199219931994199519961997199819992000200120022003200420052006200720082009201020112012201320142015201620172018201920202021202220232024202520262027202820292030203120322033203420352036203720382039204020412042204320442045204620472048204920502051205220532054205520562057205820592060206120622063206420652066206720682069207020712072207320742075207620772078207920802081208220832084208520862087208820892090209120922093209420952096209720982099210021012102210321042105210621072108210921102111211221132114211521162117211821192120212121222123212421252126212721282129213021312132213321342135213621372138213921402141214221432144214521462147214821492150215121522153215421552156215721582159216021612162216321642165216621672168216921702171217221732174217521762177217821792180218121822183218421852186218721882189219021912192219321942195219621972198219922002201220222032204220522062207220822092210221122122213221422152216221722182219222022212222222322242225222622272228222922302231223222332234223522362237223822392240224122422243224422452246224722482249225022512252225322542255225622572258225922602261226222632264226522662267226822692270227122722273227422752276227722782279228022812282228322842285228622872288228922902291229222932294229522962297229822992300230123022303230423052306230723082309231023112312231323142315231623172318231923202321232223232324232523262327232823292330233123322333233423352336233723382339234023412342234323442345234623472348234923502351235223532354235523562357235823592360236123622363236423652366236723682369237023712372237323742375237623772378237923802381238223832384238523862387238823892390239123922393239423952396239723982399240024012402240324042405240624072408240924102411241224132414241524162417241824192420242124222423242424252426242724282429243024312432243324342435243624372438243924402441244224432444244524462447244824492450245124522453245424552456245724582459246024612462246324642465246624672468246924702471247224732474247524762477247824792480248124822483248424852486248724882489249024912492249324942495249624972498249925002501250225032504250525062507250825092510251125122513251425152516251725182519252025212522252325242525252625272528252925302531253225332534253525362537253825392540254125422543254425452546254725482549255025512552255325542555255625572558255925602561256225632564256525662567256825692570257125722573257425752576257725782579258025812582258325842585258625872588258925902591259225932594259525962597259825992600260126022603260426052606260726082609261026112612261326142615261626172618261926202621262226232624262526262627262826292630263126322633263426352636263726382639264026412642264326442645264626472648264926502651265226532654265526562657265826592660266126622663266426652666266726682669267026712672267326742675267626772678267926802681268226832684268526862687268826892690269126922693269426952696269726982699270027012702270327042705270627072708270927102711271227132714271527162717271827192720272127222723272427252726272727282729273027312732273327342735273627372738273927402741274227432744274527462747274827492750275127522753275427552756275727582759276027612762276327642765276627672768276927702771277227732774277527762777277827792780278127822783278427852786278727882789279027912792279327942795279627972798279928002801280228032804280528062807280828092810281128122813281428152816281728182819282028212822282328242825282628272828282928302831283228332834283528362837283828392840284128422843284428452846284728482849285028512852285328542855285628572858285928602861286228632864286528662867286828692870287128722873287428752876287728782879288028812882288328842885288628872888288928902891289228932894289528962897289828992900290129022903290429052906290729082909291029112912291329142915291629172918291929202921292229232924292529262927292829292930293129322933293429352936293729382939294029412942294329442945294629472948294929502951295229532954295529562957295829592960296129622963296429652966296729682969297029712972297329742975297629772978297929802981298229832984298529862987298829892990299129922993299429952996299729982999300030013002300330043005300630073008300930103011301230133014301530163017301830193020302130223023302430253026302730283029303030313032303330343035303630373038303930403041304230433044304530463047304830493050305130523053305430553056305730583059306030613062
  1. ///////////////////////////////////////////////////////////////////////////////
  2. // Copyright Christopher Kormanyos 2002 - 2012.
  3. // Copyright 2012 John Maddock. Distributed under the Boost
  4. // Software License, Version 1.0. (See accompanying file
  5. // LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
  6. //
  7. // This work is based on an earlier work:
  8. // "Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations",
  9. // in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
  10. //
  11. // Note that there are no "noexcept" specifications on the functions in this file: there are too many
  12. // calls to lexical_cast (and similar) to easily analyse the code for correctness. So until compilers
  13. // can detect noexcept misuse at compile time, the only realistic option is to simply not use it here.
  14. //
  15. #ifndef BOOST_MP_CPP_DEC_FLOAT_BACKEND_HPP
  16. #define BOOST_MP_CPP_DEC_FLOAT_BACKEND_HPP
  17. #include <boost/config.hpp>
  18. #include <boost/cstdint.hpp>
  19. #include <limits>
  20. #ifndef BOOST_NO_CXX11_HDR_ARRAY
  21. #include <array>
  22. #else
  23. #include <boost/array.hpp>
  24. #endif
  25. #include <boost/cstdint.hpp>
  26. #include <boost/multiprecision/number.hpp>
  27. #include <boost/multiprecision/detail/big_lanczos.hpp>
  28. #include <boost/multiprecision/detail/dynamic_array.hpp>
  29. //
  30. // Headers required for Boost.Math integration:
  31. //
  32. #include <boost/math/policies/policy.hpp>
  33. namespace boost{
  34. namespace multiprecision{
  35. namespace backends{
  36. template <unsigned Digits10, class ExponentType = boost::int32_t, class Allocator = void>
  37. class cpp_dec_float;
  38. } // namespace
  39. template <unsigned Digits10, class ExponentType, class Allocator>
  40. struct number_category<backends::cpp_dec_float<Digits10, ExponentType, Allocator> > : public mpl::int_<number_kind_floating_point>{};
  41. namespace backends{
  42. template <unsigned Digits10, class ExponentType, class Allocator>
  43. class cpp_dec_float
  44. {
  45. private:
  46. static const boost::int32_t cpp_dec_float_digits10_setting = Digits10;
  47. // We need at least 16-bits in the exponent type to do anything sensible:
  48. BOOST_STATIC_ASSERT_MSG(boost::is_signed<ExponentType>::value, "ExponentType must be a signed built in integer type.");
  49. BOOST_STATIC_ASSERT_MSG(sizeof(ExponentType) > 1, "ExponentType is too small.");
  50. public:
  51. typedef mpl::list<long long> signed_types;
  52. typedef mpl::list<unsigned long long> unsigned_types;
  53. typedef mpl::list<long double> float_types;
  54. typedef ExponentType exponent_type;
  55. static const boost::int32_t cpp_dec_float_radix = 10L;
  56. static const boost::int32_t cpp_dec_float_digits10_limit_lo = 9L;
  57. static const boost::int32_t cpp_dec_float_digits10_limit_hi = boost::integer_traits<boost::int32_t>::const_max - 100;
  58. static const boost::int32_t cpp_dec_float_digits10 = ((cpp_dec_float_digits10_setting < cpp_dec_float_digits10_limit_lo) ? cpp_dec_float_digits10_limit_lo : ((cpp_dec_float_digits10_setting > cpp_dec_float_digits10_limit_hi) ? cpp_dec_float_digits10_limit_hi : cpp_dec_float_digits10_setting));
  59. static const ExponentType cpp_dec_float_max_exp10 = (static_cast<ExponentType>(1) << (std::numeric_limits<ExponentType>::digits - 5));
  60. static const ExponentType cpp_dec_float_min_exp10 = -cpp_dec_float_max_exp10;
  61. static const ExponentType cpp_dec_float_max_exp = static_cast<ExponentType>((cpp_dec_float_max_exp10 / 301LL) * 1000LL);
  62. static const ExponentType cpp_dec_float_min_exp = static_cast<ExponentType>((cpp_dec_float_min_exp10 / 301LL) * 1000LL);
  63. BOOST_STATIC_ASSERT((cpp_dec_float<Digits10, ExponentType, Allocator>::cpp_dec_float_max_exp10 == -cpp_dec_float<Digits10, ExponentType, Allocator>::cpp_dec_float_min_exp10));
  64. private:
  65. static const boost::int32_t cpp_dec_float_elem_digits10 = 8L;
  66. static const boost::int32_t cpp_dec_float_elem_mask = 100000000L;
  67. BOOST_STATIC_ASSERT(0 == cpp_dec_float_max_exp10 % cpp_dec_float_elem_digits10);
  68. // There are three guard limbs.
  69. // 1) The first limb has 'play' from 1...8 decimal digits.
  70. // 2) The last limb also has 'play' from 1...8 decimal digits.
  71. // 3) One limb can get lost when justifying after multiply,
  72. // as only half of the triangle is multiplied and a carry
  73. // from below is missing.
  74. static const boost::int32_t cpp_dec_float_elem_number_request = static_cast<boost::int32_t>((cpp_dec_float_digits10 / cpp_dec_float_elem_digits10) + (((cpp_dec_float_digits10 % cpp_dec_float_elem_digits10) != 0) ? 1 : 0));
  75. // The number of elements needed (with a minimum of two) plus three added guard limbs.
  76. static const boost::int32_t cpp_dec_float_elem_number = static_cast<boost::int32_t>(((cpp_dec_float_elem_number_request < 2L) ? 2L : cpp_dec_float_elem_number_request) + 3L);
  77. public:
  78. static const boost::int32_t cpp_dec_float_total_digits10 = static_cast<boost::int32_t>(cpp_dec_float_elem_number * cpp_dec_float_elem_digits10);
  79. private:
  80. typedef enum enum_fpclass_type
  81. {
  82. cpp_dec_float_finite,
  83. cpp_dec_float_inf,
  84. cpp_dec_float_NaN
  85. }
  86. fpclass_type;
  87. #ifndef BOOST_NO_CXX11_HDR_ARRAY
  88. typedef typename mpl::if_<is_void<Allocator>,
  89. std::array<boost::uint32_t, cpp_dec_float_elem_number>,
  90. detail::dynamic_array<boost::uint32_t, cpp_dec_float_elem_number, Allocator>
  91. >::type array_type;
  92. #else
  93. typedef typename mpl::if_<is_void<Allocator>,
  94. boost::array<boost::uint32_t, cpp_dec_float_elem_number>,
  95. detail::dynamic_array<boost::uint32_t, cpp_dec_float_elem_number, Allocator>
  96. >::type array_type;
  97. #endif
  98. array_type data;
  99. ExponentType exp;
  100. bool neg;
  101. fpclass_type fpclass;
  102. boost::int32_t prec_elem;
  103. //
  104. // Special values constructor:
  105. //
  106. cpp_dec_float(fpclass_type c) :
  107. data(),
  108. exp (static_cast<ExponentType>(0)),
  109. neg (false),
  110. fpclass (c),
  111. prec_elem(cpp_dec_float_elem_number) { }
  112. //
  113. // Static data initializer:
  114. //
  115. struct initializer
  116. {
  117. initializer()
  118. {
  119. cpp_dec_float<Digits10, ExponentType, Allocator>::nan();
  120. cpp_dec_float<Digits10, ExponentType, Allocator>::inf();
  121. (cpp_dec_float<Digits10, ExponentType, Allocator>::min)();
  122. (cpp_dec_float<Digits10, ExponentType, Allocator>::max)();
  123. cpp_dec_float<Digits10, ExponentType, Allocator>::zero();
  124. cpp_dec_float<Digits10, ExponentType, Allocator>::one();
  125. cpp_dec_float<Digits10, ExponentType, Allocator>::two();
  126. cpp_dec_float<Digits10, ExponentType, Allocator>::half();
  127. cpp_dec_float<Digits10, ExponentType, Allocator>::double_min();
  128. cpp_dec_float<Digits10, ExponentType, Allocator>::double_max();
  129. cpp_dec_float<Digits10, ExponentType, Allocator>::long_double_max();
  130. cpp_dec_float<Digits10, ExponentType, Allocator>::long_double_min();
  131. cpp_dec_float<Digits10, ExponentType, Allocator>::long_long_max();
  132. cpp_dec_float<Digits10, ExponentType, Allocator>::long_long_min();
  133. cpp_dec_float<Digits10, ExponentType, Allocator>::ulong_long_max();
  134. cpp_dec_float<Digits10, ExponentType, Allocator>::eps();
  135. cpp_dec_float<Digits10, ExponentType, Allocator>::pow2(0);
  136. }
  137. void do_nothing(){}
  138. };
  139. static initializer init;
  140. public:
  141. // Constructors
  142. cpp_dec_float() :
  143. data(),
  144. exp (static_cast<ExponentType>(0)),
  145. neg (false),
  146. fpclass (cpp_dec_float_finite),
  147. prec_elem(cpp_dec_float_elem_number) { }
  148. cpp_dec_float(const char* s) :
  149. data(),
  150. exp (static_cast<ExponentType>(0)),
  151. neg (false),
  152. fpclass (cpp_dec_float_finite),
  153. prec_elem(cpp_dec_float_elem_number)
  154. {
  155. *this = s;
  156. }
  157. template<class I>
  158. cpp_dec_float(I i, typename enable_if<is_unsigned<I> >::type* = 0) :
  159. data(),
  160. exp (static_cast<ExponentType>(0)),
  161. neg (false),
  162. fpclass (cpp_dec_float_finite),
  163. prec_elem(cpp_dec_float_elem_number)
  164. {
  165. from_unsigned_long_long(i);
  166. }
  167. template <class I>
  168. cpp_dec_float(I i, typename enable_if<is_signed<I> >::type* = 0) :
  169. data(),
  170. exp (static_cast<ExponentType>(0)),
  171. neg (false),
  172. fpclass (cpp_dec_float_finite),
  173. prec_elem(cpp_dec_float_elem_number)
  174. {
  175. if(i < 0)
  176. {
  177. from_unsigned_long_long(-i);
  178. negate();
  179. }
  180. else
  181. from_unsigned_long_long(i);
  182. }
  183. cpp_dec_float(const cpp_dec_float& f) :
  184. data (f.data),
  185. exp (f.exp),
  186. neg (f.neg),
  187. fpclass (f.fpclass),
  188. prec_elem(f.prec_elem) { }
  189. template <unsigned D, class ET, class A>
  190. cpp_dec_float(const cpp_dec_float<D, ET, A>& f, typename enable_if_c<D <= Digits10>::type* = 0) :
  191. data(),
  192. exp (f.exp),
  193. neg (f.neg),
  194. fpclass (static_cast<fpclass_type>(static_cast<int>(f.fpclass))),
  195. prec_elem(cpp_dec_float_elem_number)
  196. {
  197. std::copy(f.data.begin(), f.data.begin() + f.prec_elem, data.begin());
  198. }
  199. template <unsigned D, class ET, class A>
  200. explicit cpp_dec_float(const cpp_dec_float<D, ET, A>& f, typename disable_if_c<D <= Digits10>::type* = 0) :
  201. data(),
  202. exp (f.exp),
  203. neg (f.neg),
  204. fpclass (static_cast<fpclass_type>(static_cast<int>(f.fpclass))),
  205. prec_elem(cpp_dec_float_elem_number)
  206. {
  207. // TODO: this doesn't round!
  208. std::copy(f.data.begin(), f.data.begin() + prec_elem, data.begin());
  209. }
  210. template <class F>
  211. cpp_dec_float(const F val, typename enable_if<is_floating_point<F> >::type* = 0) :
  212. data(),
  213. exp (static_cast<ExponentType>(0)),
  214. neg (false),
  215. fpclass (cpp_dec_float_finite),
  216. prec_elem(cpp_dec_float_elem_number)
  217. {
  218. *this = val;
  219. }
  220. cpp_dec_float(const double val, ExponentType exponent);
  221. // Specific special values.
  222. static const cpp_dec_float& nan()
  223. {
  224. static const cpp_dec_float val(cpp_dec_float_NaN);
  225. init.do_nothing();
  226. return val;
  227. }
  228. static const cpp_dec_float& inf()
  229. {
  230. static const cpp_dec_float val(cpp_dec_float_inf);
  231. init.do_nothing();
  232. return val;
  233. }
  234. static const cpp_dec_float& (max)()
  235. {
  236. init.do_nothing();
  237. static bool init = false;
  238. static const std::string str_max = std::string("9." + std::string(static_cast<std::size_t>(cpp_dec_float_total_digits10), static_cast<char>('9')))
  239. + std::string("e+" + boost::lexical_cast<std::string>(cpp_dec_float_max_exp10));
  240. static cpp_dec_float val_max;
  241. if(!init)
  242. {
  243. init = true;
  244. val_max = str_max.c_str();
  245. }
  246. return val_max;
  247. }
  248. static const cpp_dec_float& (min)()
  249. {
  250. init.do_nothing();
  251. static bool init = false;
  252. static cpp_dec_float val_min;
  253. if(!init)
  254. {
  255. init = true;
  256. val_min = std::string("1.0e" + boost::lexical_cast<std::string>(cpp_dec_float_min_exp10)).c_str();
  257. }
  258. return val_min;
  259. }
  260. static const cpp_dec_float& zero()
  261. {
  262. init.do_nothing();
  263. static cpp_dec_float val(static_cast<unsigned long long>(0u));
  264. return val;
  265. }
  266. static const cpp_dec_float& one()
  267. {
  268. init.do_nothing();
  269. static cpp_dec_float val(static_cast<unsigned long long>(1u));
  270. return val;
  271. }
  272. static const cpp_dec_float& two()
  273. {
  274. init.do_nothing();
  275. static cpp_dec_float val(static_cast<unsigned long long>(2u));
  276. return val;
  277. }
  278. static const cpp_dec_float& half()
  279. {
  280. init.do_nothing();
  281. static cpp_dec_float val(0.5L);
  282. return val;
  283. }
  284. static const cpp_dec_float& double_min()
  285. {
  286. init.do_nothing();
  287. static cpp_dec_float val(static_cast<long double>((std::numeric_limits<double>::min)()));
  288. return val;
  289. }
  290. static const cpp_dec_float& double_max()
  291. {
  292. init.do_nothing();
  293. static cpp_dec_float val(static_cast<long double>((std::numeric_limits<double>::max)()));
  294. return val;
  295. }
  296. static const cpp_dec_float& long_double_min()
  297. {
  298. init.do_nothing();
  299. static cpp_dec_float val((std::numeric_limits<long double>::min)());
  300. return val;
  301. }
  302. static const cpp_dec_float& long_double_max()
  303. {
  304. init.do_nothing();
  305. static cpp_dec_float val((std::numeric_limits<long double>::max)());
  306. return val;
  307. }
  308. static const cpp_dec_float& long_long_max()
  309. {
  310. init.do_nothing();
  311. static cpp_dec_float val((std::numeric_limits<long long>::max)());
  312. return val;
  313. }
  314. static const cpp_dec_float& long_long_min()
  315. {
  316. init.do_nothing();
  317. static cpp_dec_float val((std::numeric_limits<long long>::min)());
  318. return val;
  319. }
  320. static const cpp_dec_float& ulong_long_max()
  321. {
  322. init.do_nothing();
  323. static cpp_dec_float val((std::numeric_limits<unsigned long long>::max)());
  324. return val;
  325. }
  326. static const cpp_dec_float& eps()
  327. {
  328. init.do_nothing();
  329. static cpp_dec_float val(1.0, 1 - (int)Digits10);
  330. return val;
  331. }
  332. // Basic operations.
  333. cpp_dec_float& operator= (const cpp_dec_float& v)
  334. {
  335. data = v.data;
  336. exp = v.exp;
  337. neg = v.neg;
  338. fpclass = v.fpclass;
  339. prec_elem = v.prec_elem;
  340. return *this;
  341. }
  342. template <unsigned D>
  343. cpp_dec_float& operator=(const cpp_dec_float<D>& f)
  344. {
  345. exp = f.exp;
  346. neg = f.neg;
  347. fpclass = static_cast<enum_fpclass_type>(static_cast<int>(f.fpclass));
  348. unsigned elems = (std::min)(f.prec_elem, cpp_dec_float_elem_number);
  349. std::copy(f.data.begin(), f.data.begin() + elems, data.begin());
  350. std::fill(data.begin() + elems, data.end(), 0);
  351. prec_elem = cpp_dec_float_elem_number;
  352. return *this;
  353. }
  354. cpp_dec_float& operator= (long long v)
  355. {
  356. if(v < 0)
  357. {
  358. from_unsigned_long_long(-v);
  359. negate();
  360. }
  361. else
  362. from_unsigned_long_long(v);
  363. return *this;
  364. }
  365. cpp_dec_float& operator= (unsigned long long v)
  366. {
  367. from_unsigned_long_long(v);
  368. return *this;
  369. }
  370. cpp_dec_float& operator= (long double v) ;
  371. cpp_dec_float& operator= (const char* v)
  372. {
  373. rd_string(v);
  374. return *this;
  375. }
  376. cpp_dec_float& operator+=(const cpp_dec_float& v) ;
  377. cpp_dec_float& operator-=(const cpp_dec_float& v) ;
  378. cpp_dec_float& operator*=(const cpp_dec_float& v) ;
  379. cpp_dec_float& operator/=(const cpp_dec_float& v) ;
  380. cpp_dec_float& add_unsigned_long_long(const unsigned long long n)
  381. {
  382. cpp_dec_float t;
  383. t.from_unsigned_long_long(n);
  384. return *this += t;
  385. }
  386. cpp_dec_float& sub_unsigned_long_long(const unsigned long long n)
  387. {
  388. cpp_dec_float t;
  389. t.from_unsigned_long_long(n);
  390. return *this -= t;
  391. }
  392. cpp_dec_float& mul_unsigned_long_long(const unsigned long long n);
  393. cpp_dec_float& div_unsigned_long_long(const unsigned long long n);
  394. // Elementary primitives.
  395. cpp_dec_float& calculate_inv (void) ;
  396. cpp_dec_float& calculate_sqrt(void) ;
  397. void negate()
  398. {
  399. if(!iszero())
  400. neg = !neg;
  401. }
  402. // Comparison functions
  403. bool isnan BOOST_PREVENT_MACRO_SUBSTITUTION(void) const { return (fpclass == cpp_dec_float_NaN); }
  404. bool isinf BOOST_PREVENT_MACRO_SUBSTITUTION(void) const { return (fpclass == cpp_dec_float_inf); }
  405. bool isfinite BOOST_PREVENT_MACRO_SUBSTITUTION(void) const { return (fpclass == cpp_dec_float_finite); }
  406. bool iszero (void) const
  407. {
  408. return ((fpclass == cpp_dec_float_finite) && (data[0u] == 0u));
  409. }
  410. bool isone (void) const ;
  411. bool isint (void) const ;
  412. bool isneg (void) const { return neg; }
  413. // Operators pre-increment and pre-decrement
  414. cpp_dec_float& operator++(void)
  415. {
  416. return *this += one();
  417. }
  418. cpp_dec_float& operator--(void)
  419. {
  420. return *this -= one();
  421. }
  422. std::string str(boost::intmax_t digits, std::ios_base::fmtflags f)const;
  423. int compare(const cpp_dec_float& v)const ;
  424. template <class V>
  425. int compare(const V& v)const
  426. {
  427. cpp_dec_float<Digits10, ExponentType, Allocator> t;
  428. t = v;
  429. return compare(t);
  430. }
  431. void swap(cpp_dec_float& v)
  432. {
  433. data.swap(v.data);
  434. std::swap(exp, v.exp);
  435. std::swap(neg, v.neg);
  436. std::swap(fpclass, v.fpclass);
  437. std::swap(prec_elem, v.prec_elem);
  438. }
  439. double extract_double (void) const;
  440. long double extract_long_double (void) const;
  441. signed long long extract_signed_long_long (void) const ;
  442. unsigned long long extract_unsigned_long_long(void) const ;
  443. void extract_parts (double& mantissa, ExponentType& exponent) const ;
  444. cpp_dec_float extract_integer_part (void) const ;
  445. void precision(const boost::int32_t prec_digits)
  446. {
  447. if(prec_digits >= cpp_dec_float_total_digits10)
  448. {
  449. prec_elem = cpp_dec_float_elem_number;
  450. }
  451. else
  452. {
  453. const boost::int32_t elems = static_cast<boost::int32_t>( static_cast<boost::int32_t>( (prec_digits + (cpp_dec_float_elem_digits10 / 2)) / cpp_dec_float_elem_digits10)
  454. + static_cast<boost::int32_t>(((prec_digits % cpp_dec_float_elem_digits10) != 0) ? 1 : 0));
  455. prec_elem = (std::min)(cpp_dec_float_elem_number, (std::max)(elems, static_cast<boost::int32_t>(2)));
  456. }
  457. }
  458. static cpp_dec_float pow2(long long i);
  459. ExponentType order()const
  460. {
  461. const bool bo_order_is_zero = ((!(isfinite)()) || (data[0] == static_cast<boost::uint32_t>(0u)));
  462. //
  463. // Binary search to find the order of the leading term:
  464. //
  465. ExponentType prefix = 0;
  466. if(data[0] >= 100000UL)
  467. {
  468. if(data[0] >= 10000000UL)
  469. {
  470. if(data[0] >= 100000000UL)
  471. {
  472. if(data[0] >= 1000000000UL)
  473. prefix = 9;
  474. else
  475. prefix = 8;
  476. }
  477. else
  478. prefix = 7;
  479. }
  480. else
  481. {
  482. if(data[0] >= 1000000UL)
  483. prefix = 6;
  484. else
  485. prefix = 5;
  486. }
  487. }
  488. else
  489. {
  490. if(data[0] >= 1000UL)
  491. {
  492. if(data[0] >= 10000UL)
  493. prefix = 4;
  494. else
  495. prefix = 3;
  496. }
  497. else
  498. {
  499. if(data[0] >= 100)
  500. prefix = 2;
  501. else if(data[0] >= 10)
  502. prefix = 1;
  503. }
  504. }
  505. return (bo_order_is_zero ? static_cast<ExponentType>(0) : static_cast<ExponentType>(exp + prefix));
  506. }
  507. template<class Archive>
  508. void serialize(Archive & ar, const unsigned int /*version*/)
  509. {
  510. for(unsigned i = 0; i < data.size(); ++i)
  511. ar & data[i];
  512. ar & exp;
  513. ar & neg;
  514. ar & fpclass;
  515. ar & prec_elem;
  516. }
  517. private:
  518. static bool data_elem_is_non_zero_predicate(const boost::uint32_t& d) { return (d != static_cast<boost::uint32_t>(0u)); }
  519. static bool data_elem_is_non_nine_predicate(const boost::uint32_t& d) { return (d != static_cast<boost::uint32_t>(cpp_dec_float::cpp_dec_float_elem_mask - 1)); }
  520. static bool char_is_nonzero_predicate(const char& c) { return (c != static_cast<char>('0')); }
  521. void from_unsigned_long_long(const unsigned long long u) ;
  522. int cmp_data(const array_type& vd) const ;
  523. static boost::uint32_t mul_loop_uv(boost::uint32_t* const u, const boost::uint32_t* const v, const boost::int32_t p) ;
  524. static boost::uint32_t mul_loop_n (boost::uint32_t* const u, boost::uint32_t n, const boost::int32_t p) ;
  525. static boost::uint32_t div_loop_n (boost::uint32_t* const u, boost::uint32_t n, const boost::int32_t p) ;
  526. bool rd_string(const char* const s);
  527. template <unsigned D, class ET, class A>
  528. friend class cpp_dec_float;
  529. };
  530. template <unsigned Digits10, class ExponentType, class Allocator>
  531. typename cpp_dec_float<Digits10, ExponentType, Allocator>::initializer cpp_dec_float<Digits10, ExponentType, Allocator>::init;
  532. template <unsigned Digits10, class ExponentType, class Allocator>
  533. const boost::int32_t cpp_dec_float<Digits10, ExponentType, Allocator>::cpp_dec_float_radix;
  534. template <unsigned Digits10, class ExponentType, class Allocator>
  535. const boost::int32_t cpp_dec_float<Digits10, ExponentType, Allocator>::cpp_dec_float_digits10_setting;
  536. template <unsigned Digits10, class ExponentType, class Allocator>
  537. const boost::int32_t cpp_dec_float<Digits10, ExponentType, Allocator>::cpp_dec_float_digits10_limit_lo;
  538. template <unsigned Digits10, class ExponentType, class Allocator>
  539. const boost::int32_t cpp_dec_float<Digits10, ExponentType, Allocator>::cpp_dec_float_digits10_limit_hi;
  540. template <unsigned Digits10, class ExponentType, class Allocator>
  541. const boost::int32_t cpp_dec_float<Digits10, ExponentType, Allocator>::cpp_dec_float_digits10;
  542. template <unsigned Digits10, class ExponentType, class Allocator>
  543. const ExponentType cpp_dec_float<Digits10, ExponentType, Allocator>::cpp_dec_float_max_exp;
  544. template <unsigned Digits10, class ExponentType, class Allocator>
  545. const ExponentType cpp_dec_float<Digits10, ExponentType, Allocator>::cpp_dec_float_min_exp;
  546. template <unsigned Digits10, class ExponentType, class Allocator>
  547. const ExponentType cpp_dec_float<Digits10, ExponentType, Allocator>::cpp_dec_float_max_exp10;
  548. template <unsigned Digits10, class ExponentType, class Allocator>
  549. const ExponentType cpp_dec_float<Digits10, ExponentType, Allocator>::cpp_dec_float_min_exp10;
  550. template <unsigned Digits10, class ExponentType, class Allocator>
  551. const boost::int32_t cpp_dec_float<Digits10, ExponentType, Allocator>::cpp_dec_float_elem_digits10;
  552. template <unsigned Digits10, class ExponentType, class Allocator>
  553. const boost::int32_t cpp_dec_float<Digits10, ExponentType, Allocator>::cpp_dec_float_elem_number_request;
  554. template <unsigned Digits10, class ExponentType, class Allocator>
  555. const boost::int32_t cpp_dec_float<Digits10, ExponentType, Allocator>::cpp_dec_float_elem_number;
  556. template <unsigned Digits10, class ExponentType, class Allocator>
  557. const boost::int32_t cpp_dec_float<Digits10, ExponentType, Allocator>::cpp_dec_float_elem_mask;
  558. template <unsigned Digits10, class ExponentType, class Allocator>
  559. cpp_dec_float<Digits10, ExponentType, Allocator>& cpp_dec_float<Digits10, ExponentType, Allocator>::operator+=(const cpp_dec_float<Digits10, ExponentType, Allocator>& v)
  560. {
  561. if((isnan)())
  562. {
  563. return *this;
  564. }
  565. if((isinf)())
  566. {
  567. if((v.isinf)() && (isneg() != v.isneg()))
  568. {
  569. *this = nan();
  570. }
  571. return *this;
  572. }
  573. if(iszero())
  574. {
  575. return operator=(v);
  576. }
  577. // Get the offset for the add/sub operation.
  578. static const ExponentType max_delta_exp = static_cast<ExponentType>((cpp_dec_float_elem_number - 1) * cpp_dec_float_elem_digits10);
  579. const ExponentType ofs_exp = static_cast<ExponentType>(exp - v.exp);
  580. // Check if the operation is out of range, requiring special handling.
  581. if(v.iszero() || (ofs_exp > max_delta_exp))
  582. {
  583. // Result is *this unchanged since v is negligible compared to *this.
  584. return *this;
  585. }
  586. else if(ofs_exp < -max_delta_exp)
  587. {
  588. // Result is *this = v since *this is negligible compared to v.
  589. return operator=(v);
  590. }
  591. // Do the add/sub operation.
  592. typename array_type::iterator p_u = data.begin();
  593. typename array_type::const_iterator p_v = v.data.begin();
  594. bool b_copy = false;
  595. const boost::int32_t ofs = static_cast<boost::int32_t>(static_cast<boost::int32_t>(ofs_exp) / cpp_dec_float_elem_digits10);
  596. array_type n_data;
  597. if(neg == v.neg)
  598. {
  599. // Add v to *this, where the data array of either *this or v
  600. // might have to be treated with a positive, negative or zero offset.
  601. // The result is stored in *this. The data are added one element
  602. // at a time, each element with carry.
  603. if(ofs >= static_cast<boost::int32_t>(0))
  604. {
  605. std::copy(v.data.begin(), v.data.end() - static_cast<size_t>(ofs), n_data.begin() + static_cast<size_t>(ofs));
  606. std::fill(n_data.begin(), n_data.begin() + static_cast<size_t>(ofs), static_cast<boost::uint32_t>(0u));
  607. p_v = n_data.begin();
  608. }
  609. else
  610. {
  611. std::copy(data.begin(), data.end() - static_cast<size_t>(-ofs), n_data.begin() + static_cast<size_t>(-ofs));
  612. std::fill(n_data.begin(), n_data.begin() + static_cast<size_t>(-ofs), static_cast<boost::uint32_t>(0u));
  613. p_u = n_data.begin();
  614. b_copy = true;
  615. }
  616. // Addition algorithm
  617. boost::uint32_t carry = static_cast<boost::uint32_t>(0u);
  618. for(boost::int32_t j = static_cast<boost::int32_t>(cpp_dec_float_elem_number - static_cast<boost::int32_t>(1)); j >= static_cast<boost::int32_t>(0); j--)
  619. {
  620. boost::uint32_t t = static_cast<boost::uint32_t>(static_cast<boost::uint32_t>(p_u[j] + p_v[j]) + carry);
  621. carry = t / static_cast<boost::uint32_t>(cpp_dec_float_elem_mask);
  622. p_u[j] = static_cast<boost::uint32_t>(t - static_cast<boost::uint32_t>(carry * static_cast<boost::uint32_t>(cpp_dec_float_elem_mask)));
  623. }
  624. if(b_copy)
  625. {
  626. data = n_data;
  627. exp = v.exp;
  628. }
  629. // There needs to be a carry into the element -1 of the array data
  630. if(carry != static_cast<boost::uint32_t>(0u))
  631. {
  632. std::copy_backward(data.begin(), data.end() - static_cast<std::size_t>(1u), data.end());
  633. data[0] = carry;
  634. exp += static_cast<ExponentType>(cpp_dec_float_elem_digits10);
  635. }
  636. }
  637. else
  638. {
  639. // Subtract v from *this, where the data array of either *this or v
  640. // might have to be treated with a positive, negative or zero offset.
  641. if((ofs > static_cast<boost::int32_t>(0))
  642. || ( (ofs == static_cast<boost::int32_t>(0))
  643. && (cmp_data(v.data) > static_cast<boost::int32_t>(0)))
  644. )
  645. {
  646. // In this case, |u| > |v| and ofs is positive.
  647. // Copy the data of v, shifted down to a lower value
  648. // into the data array m_n. Set the operand pointer p_v
  649. // to point to the copied, shifted data m_n.
  650. std::copy(v.data.begin(), v.data.end() - static_cast<size_t>(ofs), n_data.begin() + static_cast<size_t>(ofs));
  651. std::fill(n_data.begin(), n_data.begin() + static_cast<size_t>(ofs), static_cast<boost::uint32_t>(0u));
  652. p_v = n_data.begin();
  653. }
  654. else
  655. {
  656. if(ofs != static_cast<boost::int32_t>(0))
  657. {
  658. // In this case, |u| < |v| and ofs is negative.
  659. // Shift the data of u down to a lower value.
  660. std::copy_backward(data.begin(), data.end() - static_cast<size_t>(-ofs), data.end());
  661. std::fill(data.begin(), data.begin() + static_cast<size_t>(-ofs), static_cast<boost::uint32_t>(0u));
  662. }
  663. // Copy the data of v into the data array n_data.
  664. // Set the u-pointer p_u to point to m_n and the
  665. // operand pointer p_v to point to the shifted
  666. // data m_data.
  667. n_data = v.data;
  668. p_u = n_data.begin();
  669. p_v = data.begin();
  670. b_copy = true;
  671. }
  672. boost::int32_t j;
  673. // Subtraction algorithm
  674. boost::int32_t borrow = static_cast<boost::int32_t>(0);
  675. for(j = static_cast<boost::int32_t>(cpp_dec_float_elem_number - static_cast<boost::int32_t>(1)); j >= static_cast<boost::int32_t>(0); j--)
  676. {
  677. boost::int32_t t = static_cast<boost::int32_t>(static_cast<boost::int32_t>( static_cast<boost::int32_t>(p_u[j])
  678. - static_cast<boost::int32_t>(p_v[j])) - borrow);
  679. // Underflow? Borrow?
  680. if(t < static_cast<boost::int32_t>(0))
  681. {
  682. // Yes, underflow and borrow
  683. t += static_cast<boost::int32_t>(cpp_dec_float_elem_mask);
  684. borrow = static_cast<boost::int32_t>(1);
  685. }
  686. else
  687. {
  688. borrow = static_cast<boost::int32_t>(0);
  689. }
  690. p_u[j] = static_cast<boost::uint32_t>(static_cast<boost::uint32_t>(t) % static_cast<boost::uint32_t>(cpp_dec_float_elem_mask));
  691. }
  692. if(b_copy)
  693. {
  694. data = n_data;
  695. exp = v.exp;
  696. neg = v.neg;
  697. }
  698. // Is it necessary to justify the data?
  699. const typename array_type::const_iterator first_nonzero_elem = std::find_if(data.begin(), data.end(), data_elem_is_non_zero_predicate);
  700. if(first_nonzero_elem != data.begin())
  701. {
  702. if(first_nonzero_elem == data.end())
  703. {
  704. // This result of the subtraction is exactly zero.
  705. // Reset the sign and the exponent.
  706. neg = false;
  707. exp = static_cast<ExponentType>(0);
  708. }
  709. else
  710. {
  711. // Justify the data
  712. const std::size_t sj = static_cast<std::size_t>(std::distance<typename array_type::const_iterator>(data.begin(), first_nonzero_elem));
  713. std::copy(data.begin() + static_cast<std::size_t>(sj), data.end(), data.begin());
  714. std::fill(data.end() - sj, data.end(), static_cast<boost::uint32_t>(0u));
  715. exp -= static_cast<ExponentType>(sj * static_cast<std::size_t>(cpp_dec_float_elem_digits10));
  716. }
  717. }
  718. }
  719. // Check for underflow.
  720. if(iszero())
  721. {
  722. return *this = zero();
  723. }
  724. bool overflow = exp >= cpp_dec_float_max_exp10;
  725. if(exp == cpp_dec_float_max_exp10)
  726. {
  727. // Check to see if we really truly have an overflow or not...
  728. if(isneg())
  729. {
  730. cpp_dec_float t(*this);
  731. t.negate();
  732. overflow = t.compare((max)()) > 0;
  733. }
  734. else
  735. {
  736. overflow = compare((max)()) > 0;
  737. }
  738. }
  739. // Check for overflow.
  740. if(overflow)
  741. {
  742. const bool b_result_is_neg = neg;
  743. *this = inf();
  744. if(b_result_is_neg)
  745. negate();
  746. }
  747. return *this;
  748. }
  749. template <unsigned Digits10, class ExponentType, class Allocator>
  750. cpp_dec_float<Digits10, ExponentType, Allocator>& cpp_dec_float<Digits10, ExponentType, Allocator>::operator-=(const cpp_dec_float<Digits10, ExponentType, Allocator>& v)
  751. {
  752. // Use *this - v = -(-*this + v).
  753. negate();
  754. *this += v;
  755. negate();
  756. return *this;
  757. }
  758. template <unsigned Digits10, class ExponentType, class Allocator>
  759. cpp_dec_float<Digits10, ExponentType, Allocator>& cpp_dec_float<Digits10, ExponentType, Allocator>::operator*=(const cpp_dec_float<Digits10, ExponentType, Allocator>& v)
  760. {
  761. // Evaluate the sign of the result.
  762. const bool b_result_is_neg = (neg != v.neg);
  763. // Artificially set the sign of the result to be positive.
  764. neg = false;
  765. // Handle special cases like zero, inf and NaN.
  766. const bool b_u_is_inf = (isinf)();
  767. const bool b_v_is_inf = (v.isinf)();
  768. const bool b_u_is_zero = iszero();
  769. const bool b_v_is_zero = v.iszero();
  770. if( ((isnan)() || (v.isnan)())
  771. || (b_u_is_inf && b_v_is_zero)
  772. || (b_v_is_inf && b_u_is_zero)
  773. )
  774. {
  775. *this = nan();
  776. return *this;
  777. }
  778. if(b_u_is_inf || b_v_is_inf)
  779. {
  780. *this = inf();
  781. if(b_result_is_neg)
  782. negate();
  783. return *this;
  784. }
  785. if(b_u_is_zero || b_v_is_zero)
  786. {
  787. return *this = zero();
  788. }
  789. // Check for overflow or underflow.
  790. const bool u_exp_is_neg = (exp < static_cast<ExponentType>(0));
  791. const bool v_exp_is_neg = (v.exp < static_cast<ExponentType>(0));
  792. if(u_exp_is_neg == v_exp_is_neg)
  793. {
  794. // Get the unsigned base-10 exponents of *this and v and...
  795. const ExponentType u_exp = ((!u_exp_is_neg) ? exp : static_cast<ExponentType>( -exp));
  796. const ExponentType v_exp = ((!v_exp_is_neg) ? v.exp : static_cast<ExponentType>(-v.exp));
  797. // Check the range of the upcoming multiplication.
  798. const bool b_result_is_out_of_range = (v_exp >= static_cast<ExponentType>(cpp_dec_float_max_exp10 - u_exp));
  799. if(b_result_is_out_of_range)
  800. {
  801. if(u_exp_is_neg)
  802. {
  803. *this = zero();
  804. }
  805. else
  806. {
  807. *this = inf();
  808. if(b_result_is_neg)
  809. negate();
  810. }
  811. return *this;
  812. }
  813. }
  814. // Set the exponent of the result.
  815. exp += v.exp;
  816. const boost::int32_t prec_mul = (std::min)(prec_elem, v.prec_elem);
  817. const boost::uint32_t carry = mul_loop_uv(data.data(), v.data.data(), prec_mul);
  818. // Handle a potential carry.
  819. if(carry != static_cast<boost::uint32_t>(0u))
  820. {
  821. exp += cpp_dec_float_elem_digits10;
  822. // Shift the result of the multiplication one element to the right...
  823. std::copy_backward(data.begin(),
  824. data.begin() + static_cast<std::size_t>(prec_elem - static_cast<boost::int32_t>(1)),
  825. data.begin() + static_cast<std::size_t>(prec_elem));
  826. // ... And insert the carry.
  827. data.front() = carry;
  828. }
  829. // Set the sign of the result.
  830. neg = b_result_is_neg;
  831. return *this;
  832. }
  833. template <unsigned Digits10, class ExponentType, class Allocator>
  834. cpp_dec_float<Digits10, ExponentType, Allocator>& cpp_dec_float<Digits10, ExponentType, Allocator>::operator/=(const cpp_dec_float<Digits10, ExponentType, Allocator>& v)
  835. {
  836. const bool u_and_v_are_finite_and_identical = ( (isfinite)()
  837. && (fpclass == v.fpclass)
  838. && (exp == v.exp)
  839. && (cmp_data(v.data) == static_cast<boost::int32_t>(0)));
  840. if(u_and_v_are_finite_and_identical)
  841. {
  842. if(neg != v.neg)
  843. {
  844. *this = one();
  845. negate();
  846. }
  847. else
  848. *this = one();
  849. return *this;
  850. }
  851. else
  852. {
  853. if(iszero())
  854. {
  855. if((v.isnan)() || v.iszero())
  856. {
  857. return *this = v;
  858. }
  859. return *this;
  860. }
  861. cpp_dec_float t(v);
  862. t.calculate_inv();
  863. return operator*=(t);
  864. }
  865. }
  866. template <unsigned Digits10, class ExponentType, class Allocator>
  867. cpp_dec_float<Digits10, ExponentType, Allocator>& cpp_dec_float<Digits10, ExponentType, Allocator>::mul_unsigned_long_long(const unsigned long long n)
  868. {
  869. // Multiply *this with a constant unsigned long long.
  870. // Evaluate the sign of the result.
  871. const bool b_neg = neg;
  872. // Artificially set the sign of the result to be positive.
  873. neg = false;
  874. // Handle special cases like zero, inf and NaN.
  875. const bool b_u_is_inf = (isinf)();
  876. const bool b_n_is_zero = (n == static_cast<boost::int32_t>(0));
  877. if((isnan)() || (b_u_is_inf && b_n_is_zero))
  878. {
  879. return (*this = nan());
  880. }
  881. if(b_u_is_inf)
  882. {
  883. *this = inf();
  884. if(b_neg)
  885. negate();
  886. return *this;
  887. }
  888. if(iszero() || b_n_is_zero)
  889. {
  890. // Multiplication by zero.
  891. return *this = zero();
  892. }
  893. if(n >= static_cast<unsigned long long>(cpp_dec_float_elem_mask))
  894. {
  895. neg = b_neg;
  896. cpp_dec_float t;
  897. t = n;
  898. return operator*=(t);
  899. }
  900. if(n == static_cast<unsigned long long>(1u))
  901. {
  902. neg = b_neg;
  903. return *this;
  904. }
  905. // Set up the multiplication loop.
  906. const boost::uint32_t nn = static_cast<boost::uint32_t>(n);
  907. const boost::uint32_t carry = mul_loop_n(data.data(), nn, prec_elem);
  908. // Handle the carry and adjust the exponent.
  909. if(carry != static_cast<boost::uint32_t>(0u))
  910. {
  911. exp += static_cast<ExponentType>(cpp_dec_float_elem_digits10);
  912. // Shift the result of the multiplication one element to the right.
  913. std::copy_backward(data.begin(),
  914. data.begin() + static_cast<std::size_t>(prec_elem - static_cast<boost::int32_t>(1)),
  915. data.begin() + static_cast<std::size_t>(prec_elem));
  916. data.front() = static_cast<boost::uint32_t>(carry);
  917. }
  918. bool overflow = exp >= cpp_dec_float_max_exp10;
  919. if(exp == cpp_dec_float_max_exp10)
  920. {
  921. // Check to see if we really truly have an overflow or not...
  922. if(isneg())
  923. {
  924. cpp_dec_float t(*this);
  925. t.negate();
  926. overflow = t.compare((max)()) > 0;
  927. }
  928. else
  929. {
  930. overflow = compare((max)()) > 0;
  931. }
  932. }
  933. if(overflow)
  934. {
  935. *this = inf();
  936. if(b_neg)
  937. negate();
  938. return *this;
  939. }
  940. // Set the sign.
  941. neg = b_neg;
  942. return *this;
  943. }
  944. template <unsigned Digits10, class ExponentType, class Allocator>
  945. cpp_dec_float<Digits10, ExponentType, Allocator>& cpp_dec_float<Digits10, ExponentType, Allocator>::div_unsigned_long_long(const unsigned long long n)
  946. {
  947. // Divide *this by a constant unsigned long long.
  948. // Evaluate the sign of the result.
  949. const bool b_neg = neg;
  950. // Artificially set the sign of the result to be positive.
  951. neg = false;
  952. // Handle special cases like zero, inf and NaN.
  953. if((isnan)())
  954. {
  955. return *this;
  956. }
  957. if((isinf)())
  958. {
  959. *this = inf();
  960. if(b_neg)
  961. negate();
  962. return *this;
  963. }
  964. if(n == static_cast<unsigned long long>(0u))
  965. {
  966. // Divide by 0.
  967. if(iszero())
  968. {
  969. *this = nan();
  970. return *this;
  971. }
  972. else
  973. {
  974. *this = inf();
  975. if(isneg())
  976. negate();
  977. return *this;
  978. }
  979. }
  980. if(iszero())
  981. {
  982. return *this;
  983. }
  984. if(n >= static_cast<unsigned long long>(cpp_dec_float_elem_mask))
  985. {
  986. neg = b_neg;
  987. cpp_dec_float t;
  988. t = n;
  989. return operator/=(t);
  990. }
  991. const boost::uint32_t nn = static_cast<boost::uint32_t>(n);
  992. if(nn > static_cast<boost::uint32_t>(1u))
  993. {
  994. // Do the division loop.
  995. const boost::uint32_t prev = div_loop_n(data.data(), nn, prec_elem);
  996. // Determine if one leading zero is in the result data.
  997. if(data[0] == static_cast<boost::uint32_t>(0u))
  998. {
  999. // Adjust the exponent
  1000. exp -= static_cast<ExponentType>(cpp_dec_float_elem_digits10);
  1001. // Shift result of the division one element to the left.
  1002. std::copy(data.begin() + static_cast<std::size_t>(1u),
  1003. data.begin() + static_cast<std::size_t>(prec_elem - static_cast<boost::int32_t>(1)),
  1004. data.begin());
  1005. data[prec_elem - static_cast<boost::int32_t>(1)] = static_cast<boost::uint32_t>(static_cast<boost::uint64_t>(prev * static_cast<boost::uint64_t>(cpp_dec_float_elem_mask)) / nn);
  1006. }
  1007. }
  1008. // Check for underflow.
  1009. if(iszero())
  1010. {
  1011. return *this = zero();
  1012. }
  1013. // Set the sign of the result.
  1014. neg = b_neg;
  1015. return *this;
  1016. }
  1017. template <unsigned Digits10, class ExponentType, class Allocator>
  1018. cpp_dec_float<Digits10, ExponentType, Allocator>& cpp_dec_float<Digits10, ExponentType, Allocator>::calculate_inv()
  1019. {
  1020. // Compute the inverse of *this.
  1021. const bool b_neg = neg;
  1022. neg = false;
  1023. // Handle special cases like zero, inf and NaN.
  1024. if(iszero())
  1025. {
  1026. *this = inf();
  1027. if(b_neg)
  1028. negate();
  1029. return *this;
  1030. }
  1031. if((isnan)())
  1032. {
  1033. return *this;
  1034. }
  1035. if((isinf)())
  1036. {
  1037. return *this = zero();
  1038. }
  1039. if(isone())
  1040. {
  1041. if(b_neg)
  1042. negate();
  1043. return *this;
  1044. }
  1045. // Save the original *this.
  1046. cpp_dec_float<Digits10, ExponentType, Allocator> x(*this);
  1047. // Generate the initial estimate using division.
  1048. // Extract the mantissa and exponent for a "manual"
  1049. // computation of the estimate.
  1050. double dd;
  1051. ExponentType ne;
  1052. x.extract_parts(dd, ne);
  1053. // Do the inverse estimate using double precision estimates of mantissa and exponent.
  1054. operator=(cpp_dec_float<Digits10, ExponentType, Allocator>(1.0 / dd, -ne));
  1055. // Compute the inverse of *this. Quadratically convergent Newton-Raphson iteration
  1056. // is used. During the iterative steps, the precision of the calculation is limited
  1057. // to the minimum required in order to minimize the run-time.
  1058. static const boost::int32_t double_digits10_minus_a_few = std::numeric_limits<double>::digits10 - 3;
  1059. for(boost::int32_t digits = double_digits10_minus_a_few; digits <= cpp_dec_float_total_digits10; digits *= static_cast<boost::int32_t>(2))
  1060. {
  1061. // Adjust precision of the terms.
  1062. precision(static_cast<boost::int32_t>((digits + 10) * static_cast<boost::int32_t>(2)));
  1063. x.precision(static_cast<boost::int32_t>((digits + 10) * static_cast<boost::int32_t>(2)));
  1064. // Next iteration.
  1065. cpp_dec_float t(*this);
  1066. t *= x;
  1067. t -= two();
  1068. t.negate();
  1069. *this *= t;
  1070. }
  1071. neg = b_neg;
  1072. prec_elem = cpp_dec_float_elem_number;
  1073. return *this;
  1074. }
  1075. template <unsigned Digits10, class ExponentType, class Allocator>
  1076. cpp_dec_float<Digits10, ExponentType, Allocator>& cpp_dec_float<Digits10, ExponentType, Allocator>::calculate_sqrt(void)
  1077. {
  1078. // Compute the square root of *this.
  1079. if(isneg() || (!(isfinite)()))
  1080. {
  1081. *this = nan();
  1082. return *this;
  1083. }
  1084. if(iszero() || isone())
  1085. {
  1086. return *this;
  1087. }
  1088. // Save the original *this.
  1089. cpp_dec_float<Digits10, ExponentType, Allocator> x(*this);
  1090. // Generate the initial estimate using division.
  1091. // Extract the mantissa and exponent for a "manual"
  1092. // computation of the estimate.
  1093. double dd;
  1094. ExponentType ne;
  1095. extract_parts(dd, ne);
  1096. // Force the exponent to be an even multiple of two.
  1097. if((ne % static_cast<ExponentType>(2)) != static_cast<ExponentType>(0))
  1098. {
  1099. ++ne;
  1100. dd /= 10.0;
  1101. }
  1102. // Setup the iteration.
  1103. // Estimate the square root using simple manipulations.
  1104. const double sqd = std::sqrt(dd);
  1105. *this = cpp_dec_float<Digits10, ExponentType, Allocator>(sqd, static_cast<ExponentType>(ne / static_cast<ExponentType>(2)));
  1106. // Estimate 1.0 / (2.0 * x0) using simple manipulations.
  1107. cpp_dec_float<Digits10, ExponentType, Allocator> vi(0.5 / sqd, static_cast<ExponentType>(-ne / static_cast<ExponentType>(2)));
  1108. // Compute the square root of x. Coupled Newton iteration
  1109. // as described in "Pi Unleashed" is used. During the
  1110. // iterative steps, the precision of the calculation is
  1111. // limited to the minimum required in order to minimize
  1112. // the run-time.
  1113. //
  1114. // Book references:
  1115. // http://www.jjj.de/pibook/pibook.html
  1116. // http://www.amazon.com/exec/obidos/tg/detail/-/3540665722/qid=1035535482/sr=8-7/ref=sr_8_7/104-3357872-6059916?v=glance&n=507846
  1117. static const boost::uint32_t double_digits10_minus_a_few = std::numeric_limits<double>::digits10 - 3;
  1118. for(boost::int32_t digits = double_digits10_minus_a_few; digits <= cpp_dec_float_total_digits10; digits *= 2u)
  1119. {
  1120. // Adjust precision of the terms.
  1121. precision((digits + 10) * 2);
  1122. vi.precision((digits + 10) * 2);
  1123. // Next iteration of vi
  1124. cpp_dec_float t(*this);
  1125. t *= vi;
  1126. t.negate();
  1127. t.mul_unsigned_long_long(2u);
  1128. t += one();
  1129. t *= vi;
  1130. vi += t;
  1131. // Next iteration of *this
  1132. t = *this;
  1133. t *= *this;
  1134. t.negate();
  1135. t += x;
  1136. t *= vi;
  1137. *this += t;
  1138. }
  1139. prec_elem = cpp_dec_float_elem_number;
  1140. return *this;
  1141. }
  1142. template <unsigned Digits10, class ExponentType, class Allocator>
  1143. int cpp_dec_float<Digits10, ExponentType, Allocator>::cmp_data(const array_type& vd) const
  1144. {
  1145. // Compare the data of *this with those of v.
  1146. // Return +1 for *this > v
  1147. // 0 for *this = v
  1148. // -1 for *this < v
  1149. const std::pair<typename array_type::const_iterator, typename array_type::const_iterator> mismatch_pair = std::mismatch(data.begin(), data.end(), vd.begin());
  1150. const bool is_equal = ((mismatch_pair.first == data.end()) && (mismatch_pair.second == vd.end()));
  1151. if(is_equal)
  1152. {
  1153. return 0;
  1154. }
  1155. else
  1156. {
  1157. return ((*mismatch_pair.first > *mismatch_pair.second) ? 1 : -1);
  1158. }
  1159. }
  1160. template <unsigned Digits10, class ExponentType, class Allocator>
  1161. int cpp_dec_float<Digits10, ExponentType, Allocator>::compare(const cpp_dec_float& v) const
  1162. {
  1163. // Compare v with *this.
  1164. // Return +1 for *this > v
  1165. // 0 for *this = v
  1166. // -1 for *this < v
  1167. // Handle all non-finite cases.
  1168. if((!(isfinite)()) || (!(v.isfinite)()))
  1169. {
  1170. // NaN can never equal NaN. Return an implementation-dependent
  1171. // signed result. Also note that comparison of NaN with NaN
  1172. // using operators greater-than or less-than is undefined.
  1173. if((isnan)() || (v.isnan)()) { return ((isnan)() ? 1 : -1); }
  1174. if((isinf)() && (v.isinf)())
  1175. {
  1176. // Both *this and v are infinite. They are equal if they have the same sign.
  1177. // Otherwise, *this is less than v if and only if *this is negative.
  1178. return ((neg == v.neg) ? 0 : (neg ? -1 : 1));
  1179. }
  1180. if((isinf)())
  1181. {
  1182. // *this is infinite, but v is finite.
  1183. // So negative infinite *this is less than any finite v.
  1184. // Whereas positive infinite *this is greater than any finite v.
  1185. return (isneg() ? -1 : 1);
  1186. }
  1187. else
  1188. {
  1189. // *this is finite, and v is infinite.
  1190. // So any finite *this is greater than negative infinite v.
  1191. // Whereas any finite *this is less than positive infinite v.
  1192. return (v.neg ? 1 : -1);
  1193. }
  1194. }
  1195. // And now handle all *finite* cases.
  1196. if(iszero())
  1197. {
  1198. // The value of *this is zero and v is either zero or non-zero.
  1199. return (v.iszero() ? 0
  1200. : (v.neg ? 1 : -1));
  1201. }
  1202. else if(v.iszero())
  1203. {
  1204. // The value of v is zero and *this is non-zero.
  1205. return (neg ? -1 : 1);
  1206. }
  1207. else
  1208. {
  1209. // Both *this and v are non-zero.
  1210. if(neg != v.neg)
  1211. {
  1212. // The signs are different.
  1213. return (neg ? -1 : 1);
  1214. }
  1215. else if(exp != v.exp)
  1216. {
  1217. // The signs are the same and the exponents are different.
  1218. const int val_cexpression = ((exp < v.exp) ? 1 : -1);
  1219. return (neg ? val_cexpression : -val_cexpression);
  1220. }
  1221. else
  1222. {
  1223. // The signs are the same and the exponents are the same.
  1224. // Compare the data.
  1225. const int val_cmp_data = cmp_data(v.data);
  1226. return ((!neg) ? val_cmp_data : -val_cmp_data);
  1227. }
  1228. }
  1229. }
  1230. template <unsigned Digits10, class ExponentType, class Allocator>
  1231. bool cpp_dec_float<Digits10, ExponentType, Allocator>::isone() const
  1232. {
  1233. // Check if the value of *this is identically 1 or very close to 1.
  1234. const bool not_negative_and_is_finite = ((!neg) && (isfinite)());
  1235. if(not_negative_and_is_finite)
  1236. {
  1237. if((data[0u] == static_cast<boost::uint32_t>(1u)) && (exp == static_cast<ExponentType>(0)))
  1238. {
  1239. const typename array_type::const_iterator it_non_zero = std::find_if(data.begin(), data.end(), data_elem_is_non_zero_predicate);
  1240. return (it_non_zero == data.end());
  1241. }
  1242. else if((data[0u] == static_cast<boost::uint32_t>(cpp_dec_float_elem_mask - 1)) && (exp == static_cast<ExponentType>(-cpp_dec_float_elem_digits10)))
  1243. {
  1244. const typename array_type::const_iterator it_non_nine = std::find_if(data.begin(), data.end(), data_elem_is_non_nine_predicate);
  1245. return (it_non_nine == data.end());
  1246. }
  1247. }
  1248. return false;
  1249. }
  1250. template <unsigned Digits10, class ExponentType, class Allocator>
  1251. bool cpp_dec_float<Digits10, ExponentType, Allocator>::isint() const
  1252. {
  1253. if(fpclass != cpp_dec_float_finite) { return false; }
  1254. if(iszero()) { return true; }
  1255. if(exp < static_cast<ExponentType>(0)) { return false; } // |*this| < 1.
  1256. const typename array_type::size_type offset_decimal_part = static_cast<typename array_type::size_type>(exp / cpp_dec_float_elem_digits10) + 1u;
  1257. if(offset_decimal_part >= static_cast<typename array_type::size_type>(cpp_dec_float_elem_number))
  1258. {
  1259. // The number is too large to resolve the integer part.
  1260. // It considered to be a pure integer.
  1261. return true;
  1262. }
  1263. typename array_type::const_iterator it_non_zero = std::find_if(data.begin() + offset_decimal_part, data.end(), data_elem_is_non_zero_predicate);
  1264. return (it_non_zero == data.end());
  1265. }
  1266. template <unsigned Digits10, class ExponentType, class Allocator>
  1267. void cpp_dec_float<Digits10, ExponentType, Allocator>::extract_parts(double& mantissa, ExponentType& exponent) const
  1268. {
  1269. // Extract the approximate parts mantissa and base-10 exponent from the input cpp_dec_float<Digits10, ExponentType, Allocator> value x.
  1270. // Extracts the mantissa and exponent.
  1271. exponent = exp;
  1272. boost::uint32_t p10 = static_cast<boost::uint32_t>(1u);
  1273. boost::uint32_t test = data[0u];
  1274. for(;;)
  1275. {
  1276. test /= static_cast<boost::uint32_t>(10u);
  1277. if(test == static_cast<boost::uint32_t>(0u))
  1278. {
  1279. break;
  1280. }
  1281. p10 *= static_cast<boost::uint32_t>(10u);
  1282. ++exponent;
  1283. }
  1284. // Establish the upper bound of limbs for extracting the double.
  1285. const int max_elem_in_double_count = static_cast<int>(static_cast<boost::int32_t>(std::numeric_limits<double>::digits10) / cpp_dec_float_elem_digits10)
  1286. + (static_cast<int>(static_cast<boost::int32_t>(std::numeric_limits<double>::digits10) % cpp_dec_float_elem_digits10) != 0 ? 1 : 0)
  1287. + 1;
  1288. // And make sure this upper bound stays within bounds of the elems.
  1289. const std::size_t max_elem_extract_count = static_cast<std::size_t>((std::min)(static_cast<boost::int32_t>(max_elem_in_double_count), cpp_dec_float_elem_number));
  1290. // Extract into the mantissa the first limb, extracted as a double.
  1291. mantissa = static_cast<double>(data[0]);
  1292. double scale = 1.0;
  1293. // Extract the rest of the mantissa piecewise from the limbs.
  1294. for(std::size_t i = 1u; i < max_elem_extract_count; i++)
  1295. {
  1296. scale /= static_cast<double>(cpp_dec_float_elem_mask);
  1297. mantissa += (static_cast<double>(data[i]) * scale);
  1298. }
  1299. mantissa /= static_cast<double>(p10);
  1300. if(neg) { mantissa = -mantissa; }
  1301. }
  1302. template <unsigned Digits10, class ExponentType, class Allocator>
  1303. double cpp_dec_float<Digits10, ExponentType, Allocator>::extract_double(void) const
  1304. {
  1305. // Returns the double conversion of a cpp_dec_float<Digits10, ExponentType, Allocator>.
  1306. // Check for non-normal cpp_dec_float<Digits10, ExponentType, Allocator>.
  1307. if(!(isfinite)())
  1308. {
  1309. if((isnan)())
  1310. {
  1311. return std::numeric_limits<double>::quiet_NaN();
  1312. }
  1313. else
  1314. {
  1315. return ((!neg) ? std::numeric_limits<double>::infinity()
  1316. : -std::numeric_limits<double>::infinity());
  1317. }
  1318. }
  1319. cpp_dec_float<Digits10, ExponentType, Allocator> xx(*this);
  1320. if(xx.isneg())
  1321. xx.negate();
  1322. // Check if *this cpp_dec_float<Digits10, ExponentType, Allocator> is zero.
  1323. if(iszero() || (xx.compare(double_min()) < 0))
  1324. {
  1325. return 0.0;
  1326. }
  1327. // Check if *this cpp_dec_float<Digits10, ExponentType, Allocator> exceeds the maximum of double.
  1328. if(xx.compare(double_max()) > 0)
  1329. {
  1330. return ((!neg) ? std::numeric_limits<double>::infinity()
  1331. : -std::numeric_limits<double>::infinity());
  1332. }
  1333. std::stringstream ss;
  1334. ss << str(std::numeric_limits<double>::digits10 + (2 + 1), std::ios_base::scientific);
  1335. double d;
  1336. ss >> d;
  1337. return d;
  1338. }
  1339. template <unsigned Digits10, class ExponentType, class Allocator>
  1340. long double cpp_dec_float<Digits10, ExponentType, Allocator>::extract_long_double(void) const
  1341. {
  1342. // Returns the long double conversion of a cpp_dec_float<Digits10, ExponentType, Allocator>.
  1343. // Check if *this cpp_dec_float<Digits10, ExponentType, Allocator> is subnormal.
  1344. if(!(isfinite)())
  1345. {
  1346. if((isnan)())
  1347. {
  1348. return std::numeric_limits<long double>::quiet_NaN();
  1349. }
  1350. else
  1351. {
  1352. return ((!neg) ? std::numeric_limits<long double>::infinity()
  1353. : -std::numeric_limits<long double>::infinity());
  1354. }
  1355. }
  1356. cpp_dec_float<Digits10, ExponentType, Allocator> xx(*this);
  1357. if(xx.isneg())
  1358. xx.negate();
  1359. // Check if *this cpp_dec_float<Digits10, ExponentType, Allocator> is zero.
  1360. if(iszero() || (xx.compare(long_double_min()) < 0))
  1361. {
  1362. return static_cast<long double>(0.0);
  1363. }
  1364. // Check if *this cpp_dec_float<Digits10, ExponentType, Allocator> exceeds the maximum of double.
  1365. if(xx.compare(long_double_max()) > 0)
  1366. {
  1367. return ((!neg) ? std::numeric_limits<long double>::infinity()
  1368. : -std::numeric_limits<long double>::infinity());
  1369. }
  1370. std::stringstream ss;
  1371. ss << str(std::numeric_limits<long double>::digits10 + (2 + 1), std::ios_base::scientific);
  1372. long double ld;
  1373. ss >> ld;
  1374. return ld;
  1375. }
  1376. template <unsigned Digits10, class ExponentType, class Allocator>
  1377. signed long long cpp_dec_float<Digits10, ExponentType, Allocator>::extract_signed_long_long(void) const
  1378. {
  1379. // Extracts a signed long long from *this.
  1380. // If (x > maximum of signed long long) or (x < minimum of signed long long),
  1381. // then the maximum or minimum of signed long long is returned accordingly.
  1382. if(exp < static_cast<ExponentType>(0))
  1383. {
  1384. return static_cast<signed long long>(0);
  1385. }
  1386. const bool b_neg = isneg();
  1387. unsigned long long val;
  1388. if((!b_neg) && (compare(long_long_max()) > 0))
  1389. {
  1390. return (std::numeric_limits<signed long long>::max)();
  1391. }
  1392. else if(b_neg && (compare(long_long_min()) < 0))
  1393. {
  1394. return (std::numeric_limits<signed long long>::min)();
  1395. }
  1396. else
  1397. {
  1398. // Extract the data into an unsigned long long value.
  1399. cpp_dec_float<Digits10, ExponentType, Allocator> xn(extract_integer_part());
  1400. if(xn.isneg())
  1401. xn.negate();
  1402. val = static_cast<unsigned long long>(xn.data[0]);
  1403. const boost::int32_t imax = (std::min)(static_cast<boost::int32_t>(static_cast<boost::int32_t>(xn.exp) / cpp_dec_float_elem_digits10), static_cast<boost::int32_t>(cpp_dec_float_elem_number - static_cast<boost::int32_t>(1)));
  1404. for(boost::int32_t i = static_cast<boost::int32_t>(1); i <= imax; i++)
  1405. {
  1406. val *= static_cast<unsigned long long>(cpp_dec_float_elem_mask);
  1407. val += static_cast<unsigned long long>(xn.data[i]);
  1408. }
  1409. }
  1410. return ((!b_neg) ? static_cast<signed long long>(val) : static_cast<signed long long>(-static_cast<signed long long>(val)));
  1411. }
  1412. template <unsigned Digits10, class ExponentType, class Allocator>
  1413. unsigned long long cpp_dec_float<Digits10, ExponentType, Allocator>::extract_unsigned_long_long(void) const
  1414. {
  1415. // Extracts an unsigned long long from *this.
  1416. // If x exceeds the maximum of unsigned long long,
  1417. // then the maximum of unsigned long long is returned.
  1418. // If x is negative, then the unsigned long long cast of
  1419. // the signed long long extracted value is returned.
  1420. if(isneg())
  1421. {
  1422. return static_cast<unsigned long long>(extract_signed_long_long());
  1423. }
  1424. if(exp < static_cast<ExponentType>(0))
  1425. {
  1426. return static_cast<unsigned long long>(0u);
  1427. }
  1428. const cpp_dec_float<Digits10, ExponentType, Allocator> xn(extract_integer_part());
  1429. unsigned long long val;
  1430. if(xn.compare(ulong_long_max()) > 0)
  1431. {
  1432. return (std::numeric_limits<unsigned long long>::max)();
  1433. }
  1434. else
  1435. {
  1436. // Extract the data into an unsigned long long value.
  1437. val = static_cast<unsigned long long>(xn.data[0]);
  1438. const boost::int32_t imax = (std::min)(static_cast<boost::int32_t>(static_cast<boost::int32_t>(xn.exp) / cpp_dec_float_elem_digits10), static_cast<boost::int32_t>(cpp_dec_float_elem_number - static_cast<boost::int32_t>(1)));
  1439. for(boost::int32_t i = static_cast<boost::int32_t>(1); i <= imax; i++)
  1440. {
  1441. val *= static_cast<unsigned long long>(cpp_dec_float_elem_mask);
  1442. val += static_cast<unsigned long long>(xn.data[i]);
  1443. }
  1444. }
  1445. return val;
  1446. }
  1447. template <unsigned Digits10, class ExponentType, class Allocator>
  1448. cpp_dec_float<Digits10, ExponentType, Allocator> cpp_dec_float<Digits10, ExponentType, Allocator>::extract_integer_part(void) const
  1449. {
  1450. // Compute the signed integer part of x.
  1451. if(!(isfinite)())
  1452. {
  1453. return *this;
  1454. }
  1455. if(exp < static_cast<ExponentType>(0))
  1456. {
  1457. // The absolute value of the number is smaller than 1.
  1458. // Thus the integer part is zero.
  1459. return zero();
  1460. }
  1461. // Truncate the digits from the decimal part, including guard digits
  1462. // that do not belong to the integer part.
  1463. // Make a local copy.
  1464. cpp_dec_float<Digits10, ExponentType, Allocator> x = *this;
  1465. // Clear out the decimal portion
  1466. const size_t first_clear = (static_cast<size_t>(x.exp) / static_cast<size_t>(cpp_dec_float_elem_digits10)) + 1u;
  1467. const size_t last_clear = static_cast<size_t>(cpp_dec_float_elem_number);
  1468. if(first_clear < last_clear)
  1469. std::fill(x.data.begin() + first_clear, x.data.begin() + last_clear, static_cast<boost::uint32_t>(0u));
  1470. return x;
  1471. }
  1472. template <unsigned Digits10, class ExponentType, class Allocator>
  1473. std::string cpp_dec_float<Digits10, ExponentType, Allocator>::str(boost::intmax_t number_of_digits, std::ios_base::fmtflags f) const
  1474. {
  1475. if((this->isinf)())
  1476. {
  1477. if(this->isneg())
  1478. return "-inf";
  1479. else if(f & std::ios_base::showpos)
  1480. return "+inf";
  1481. else
  1482. return "inf";
  1483. }
  1484. else if((this->isnan)())
  1485. {
  1486. return "nan";
  1487. }
  1488. std::string str;
  1489. boost::intmax_t org_digits(number_of_digits);
  1490. ExponentType my_exp = order();
  1491. if(number_of_digits == 0)
  1492. number_of_digits = cpp_dec_float_total_digits10;
  1493. if(f & std::ios_base::fixed)
  1494. {
  1495. number_of_digits += my_exp + 1;
  1496. }
  1497. else if(f & std::ios_base::scientific)
  1498. ++number_of_digits;
  1499. // Determine the number of elements needed to provide the requested digits from cpp_dec_float<Digits10, ExponentType, Allocator>.
  1500. const std::size_t number_of_elements = (std::min)(static_cast<std::size_t>((number_of_digits / static_cast<std::size_t>(cpp_dec_float_elem_digits10)) + 2u),
  1501. static_cast<std::size_t>(cpp_dec_float_elem_number));
  1502. // Extract the remaining digits from cpp_dec_float<Digits10, ExponentType, Allocator> after the decimal point.
  1503. str = boost::lexical_cast<std::string>(data[0]);
  1504. // Extract all of the digits from cpp_dec_float<Digits10, ExponentType, Allocator>, beginning with the first data element.
  1505. for(std::size_t i = static_cast<std::size_t>(1u); i < number_of_elements; i++)
  1506. {
  1507. std::stringstream ss;
  1508. ss << std::setw(static_cast<std::streamsize>(cpp_dec_float_elem_digits10))
  1509. << std::setfill(static_cast<char>('0'))
  1510. << data[i];
  1511. str += ss.str();
  1512. }
  1513. bool have_leading_zeros = false;
  1514. if(number_of_digits == 0)
  1515. {
  1516. // We only get here if the output format is "fixed" and we just need to
  1517. // round the first non-zero digit.
  1518. number_of_digits -= my_exp + 1; // reset to original value
  1519. str.insert(0, std::string::size_type(number_of_digits), '0');
  1520. have_leading_zeros = true;
  1521. }
  1522. if(number_of_digits < 0)
  1523. {
  1524. str = "0";
  1525. if(isneg())
  1526. str.insert(0, 1, '-');
  1527. boost::multiprecision::detail::format_float_string(str, 0, number_of_digits - my_exp - 1, f, this->iszero());
  1528. return str;
  1529. }
  1530. else
  1531. {
  1532. // Cut the output to the size of the precision.
  1533. if(str.length() > static_cast<std::string::size_type>(number_of_digits))
  1534. {
  1535. // Get the digit after the last needed digit for rounding
  1536. const boost::uint32_t round = static_cast<boost::uint32_t>(static_cast<boost::uint32_t>(str[static_cast<std::string::size_type>(number_of_digits)]) - static_cast<boost::uint32_t>('0'));
  1537. bool need_round_up = round >= 5u;
  1538. if(round == 5u)
  1539. {
  1540. const boost::uint32_t ix = static_cast<boost::uint32_t>(static_cast<boost::uint32_t>(str[static_cast<std::string::size_type>(number_of_digits - 1)]) - static_cast<boost::uint32_t>('0'));
  1541. if((ix & 1u) == 0)
  1542. {
  1543. // We have an even digit followed by a 5, so we might not actually need to round up
  1544. // if all the remaining digits are zero:
  1545. if(str.find_first_not_of('0', static_cast<std::string::size_type>(number_of_digits + 1)) == std::string::npos)
  1546. {
  1547. bool all_zeros = true;
  1548. // No none-zero trailing digits in the string, now check whatever parts we didn't convert to the string:
  1549. for(std::size_t i = number_of_elements; i < data.size(); i++)
  1550. {
  1551. if(data[i])
  1552. {
  1553. all_zeros = false;
  1554. break;
  1555. }
  1556. }
  1557. if(all_zeros)
  1558. need_round_up = false; // tie break - round to even.
  1559. }
  1560. }
  1561. }
  1562. // Truncate the string
  1563. str.erase(static_cast<std::string::size_type>(number_of_digits));
  1564. if(need_round_up)
  1565. {
  1566. std::size_t ix = static_cast<std::size_t>(str.length() - 1u);
  1567. // Every trailing 9 must be rounded up
  1568. while(ix && (static_cast<boost::int32_t>(str.at(ix)) - static_cast<boost::int32_t>('0') == static_cast<boost::int32_t>(9)))
  1569. {
  1570. str.at(ix) = static_cast<char>('0');
  1571. --ix;
  1572. }
  1573. if(!ix)
  1574. {
  1575. // There were nothing but trailing nines.
  1576. if(static_cast<boost::int32_t>(static_cast<boost::int32_t>(str.at(ix)) - static_cast<boost::int32_t>(0x30)) == static_cast<boost::int32_t>(9))
  1577. {
  1578. // Increment up to the next order and adjust exponent.
  1579. str.at(ix) = static_cast<char>('1');
  1580. ++my_exp;
  1581. }
  1582. else
  1583. {
  1584. // Round up this digit.
  1585. ++str.at(ix);
  1586. }
  1587. }
  1588. else
  1589. {
  1590. // Round up the last digit.
  1591. ++str[ix];
  1592. }
  1593. }
  1594. }
  1595. }
  1596. if(have_leading_zeros)
  1597. {
  1598. // We need to take the zeros back out again, and correct the exponent
  1599. // if we rounded up:
  1600. if(str[std::string::size_type(number_of_digits - 1)] != '0')
  1601. {
  1602. ++my_exp;
  1603. str.erase(0, std::string::size_type(number_of_digits - 1));
  1604. }
  1605. else
  1606. str.erase(0, std::string::size_type(number_of_digits));
  1607. }
  1608. if(isneg())
  1609. str.insert(0, 1, '-');
  1610. boost::multiprecision::detail::format_float_string(str, my_exp, org_digits, f, this->iszero());
  1611. return str;
  1612. }
  1613. template <unsigned Digits10, class ExponentType, class Allocator>
  1614. bool cpp_dec_float<Digits10, ExponentType, Allocator>::rd_string(const char* const s)
  1615. {
  1616. try{
  1617. std::string str(s);
  1618. // TBD: Using several regular expressions may significantly reduce
  1619. // the code complexity (and perhaps the run-time) of rd_string().
  1620. // Get a possible exponent and remove it.
  1621. exp = static_cast<ExponentType>(0);
  1622. std::size_t pos;
  1623. if( ((pos = str.find('e')) != std::string::npos)
  1624. || ((pos = str.find('E')) != std::string::npos)
  1625. )
  1626. {
  1627. // Remove the exponent part from the string.
  1628. exp = boost::lexical_cast<ExponentType>(static_cast<const char*>(str.c_str() + (pos + 1u)));
  1629. str = str.substr(static_cast<std::size_t>(0u), pos);
  1630. }
  1631. // Get a possible +/- sign and remove it.
  1632. neg = false;
  1633. if(str.size())
  1634. {
  1635. if(str[0] == '-')
  1636. {
  1637. neg = true;
  1638. str.erase(0, 1);
  1639. }
  1640. else if(str[0] == '+')
  1641. {
  1642. str.erase(0, 1);
  1643. }
  1644. }
  1645. //
  1646. // Special cases for infinities and NaN's:
  1647. //
  1648. if((str == "inf") || (str == "INF") || (str == "infinity") || (str == "INFINITY"))
  1649. {
  1650. if(neg)
  1651. {
  1652. *this = this->inf();
  1653. this->negate();
  1654. }
  1655. else
  1656. *this = this->inf();
  1657. return true;
  1658. }
  1659. if((str.size() >= 3) && ((str.substr(0, 3) == "nan") || (str.substr(0, 3) == "NAN") || (str.substr(0, 3) == "NaN")))
  1660. {
  1661. *this = this->nan();
  1662. return true;
  1663. }
  1664. // Remove the leading zeros for all input types.
  1665. const std::string::iterator fwd_it_leading_zero = std::find_if(str.begin(), str.end(), char_is_nonzero_predicate);
  1666. if(fwd_it_leading_zero != str.begin())
  1667. {
  1668. if(fwd_it_leading_zero == str.end())
  1669. {
  1670. // The string contains nothing but leading zeros.
  1671. // This string represents zero.
  1672. operator=(zero());
  1673. return true;
  1674. }
  1675. else
  1676. {
  1677. str.erase(str.begin(), fwd_it_leading_zero);
  1678. }
  1679. }
  1680. // Put the input string into the standard cpp_dec_float<Digits10, ExponentType, Allocator> input form
  1681. // aaa.bbbbE+/-n, where aaa has 1...cpp_dec_float_elem_digits10, bbbb has an
  1682. // even multiple of cpp_dec_float_elem_digits10 which are possibly zero padded
  1683. // on the right-end, and n is a signed 64-bit integer which is an
  1684. // even multiple of cpp_dec_float_elem_digits10.
  1685. // Find a possible decimal point.
  1686. pos = str.find(static_cast<char>('.'));
  1687. if(pos != std::string::npos)
  1688. {
  1689. // Remove all trailing insignificant zeros.
  1690. const std::string::const_reverse_iterator rit_non_zero = std::find_if(str.rbegin(), str.rend(), char_is_nonzero_predicate);
  1691. if(rit_non_zero != str.rbegin())
  1692. {
  1693. const std::string::size_type ofs = str.length() - std::distance<std::string::const_reverse_iterator>(str.rbegin(), rit_non_zero);
  1694. str.erase(str.begin() + ofs, str.end());
  1695. }
  1696. // Check if the input is identically zero.
  1697. if(str == std::string("."))
  1698. {
  1699. operator=(zero());
  1700. return true;
  1701. }
  1702. // Remove leading significant zeros just after the decimal point
  1703. // and adjust the exponent accordingly.
  1704. // Note that the while-loop operates only on strings of the form ".000abcd..."
  1705. // and peels away the zeros just after the decimal point.
  1706. if(str.at(static_cast<std::size_t>(0u)) == static_cast<char>('.'))
  1707. {
  1708. const std::string::iterator it_non_zero = std::find_if(str.begin() + 1u, str.end(), char_is_nonzero_predicate);
  1709. std::size_t delta_exp = static_cast<std::size_t>(0u);
  1710. if(str.at(static_cast<std::size_t>(1u)) == static_cast<char>('0'))
  1711. {
  1712. delta_exp = std::distance<std::string::const_iterator>(str.begin() + 1u, it_non_zero);
  1713. }
  1714. // Bring one single digit into the mantissa and adjust the exponent accordingly.
  1715. str.erase(str.begin(), it_non_zero);
  1716. str.insert(static_cast<std::size_t>(1u), ".");
  1717. exp -= static_cast<ExponentType>(delta_exp + 1u);
  1718. }
  1719. }
  1720. else
  1721. {
  1722. // Input string has no decimal point: Append decimal point.
  1723. str.append(".");
  1724. }
  1725. // Shift the decimal point such that the exponent is an even multiple of cpp_dec_float_elem_digits10.
  1726. std::size_t n_shift = static_cast<std::size_t>(0u);
  1727. const std::size_t n_exp_rem = static_cast<std::size_t>(exp % static_cast<ExponentType>(cpp_dec_float_elem_digits10));
  1728. if((exp % static_cast<ExponentType>(cpp_dec_float_elem_digits10)) != static_cast<ExponentType>(0))
  1729. {
  1730. n_shift = ((exp < static_cast<ExponentType>(0))
  1731. ? static_cast<std::size_t>(n_exp_rem + static_cast<std::size_t>(cpp_dec_float_elem_digits10))
  1732. : static_cast<std::size_t>(n_exp_rem));
  1733. }
  1734. // Make sure that there are enough digits for the decimal point shift.
  1735. pos = str.find(static_cast<char>('.'));
  1736. std::size_t pos_plus_one = static_cast<std::size_t>(pos + 1u);
  1737. if((str.length() - pos_plus_one) < n_shift)
  1738. {
  1739. const std::size_t sz = static_cast<std::size_t>(n_shift - (str.length() - pos_plus_one));
  1740. str.append(std::string(sz, static_cast<char>('0')));
  1741. }
  1742. // Do the decimal point shift.
  1743. if(n_shift != static_cast<std::size_t>(0u))
  1744. {
  1745. str.insert(static_cast<std::size_t>(pos_plus_one + n_shift), ".");
  1746. str.erase(pos, static_cast<std::size_t>(1u));
  1747. exp -= static_cast<ExponentType>(n_shift);
  1748. }
  1749. // Cut the size of the mantissa to <= cpp_dec_float_elem_digits10.
  1750. pos = str.find(static_cast<char>('.'));
  1751. pos_plus_one = static_cast<std::size_t>(pos + 1u);
  1752. if(pos > static_cast<std::size_t>(cpp_dec_float_elem_digits10))
  1753. {
  1754. const boost::int32_t n_pos = static_cast<boost::int32_t>(pos);
  1755. const boost::int32_t n_rem_is_zero = ((static_cast<boost::int32_t>(n_pos % cpp_dec_float_elem_digits10) == static_cast<boost::int32_t>(0)) ? static_cast<boost::int32_t>(1) : static_cast<boost::int32_t>(0));
  1756. const boost::int32_t n = static_cast<boost::int32_t>(static_cast<boost::int32_t>(n_pos / cpp_dec_float_elem_digits10) - n_rem_is_zero);
  1757. str.insert(static_cast<std::size_t>(static_cast<boost::int32_t>(n_pos - static_cast<boost::int32_t>(n * cpp_dec_float_elem_digits10))), ".");
  1758. str.erase(pos_plus_one, static_cast<std::size_t>(1u));
  1759. exp += static_cast<ExponentType>(static_cast<ExponentType>(n) * static_cast<ExponentType>(cpp_dec_float_elem_digits10));
  1760. }
  1761. // Pad the decimal part such that its value is an even
  1762. // multiple of cpp_dec_float_elem_digits10.
  1763. pos = str.find(static_cast<char>('.'));
  1764. pos_plus_one = static_cast<std::size_t>(pos + 1u);
  1765. const boost::int32_t n_dec = static_cast<boost::int32_t>(static_cast<boost::int32_t>(str.length() - 1u) - static_cast<boost::int32_t>(pos));
  1766. const boost::int32_t n_rem = static_cast<boost::int32_t>(n_dec % cpp_dec_float_elem_digits10);
  1767. boost::int32_t n_cnt = ((n_rem != static_cast<boost::int32_t>(0)) ? static_cast<boost::int32_t>(cpp_dec_float_elem_digits10 - n_rem)
  1768. : static_cast<boost::int32_t>(0));
  1769. if(n_cnt != static_cast<boost::int32_t>(0))
  1770. {
  1771. str.append(static_cast<std::size_t>(n_cnt), static_cast<char>('0'));
  1772. }
  1773. // Truncate decimal part if it is too long.
  1774. const std::size_t max_dec = static_cast<std::size_t>((cpp_dec_float_elem_number - 1) * cpp_dec_float_elem_digits10);
  1775. if(static_cast<std::size_t>(str.length() - pos) > max_dec)
  1776. {
  1777. str = str.substr(static_cast<std::size_t>(0u),
  1778. static_cast<std::size_t>(pos_plus_one + max_dec));
  1779. }
  1780. // Now the input string has the standard cpp_dec_float<Digits10, ExponentType, Allocator> input form.
  1781. // (See the comment above.)
  1782. // Set all the data elements to 0.
  1783. std::fill(data.begin(), data.end(), static_cast<boost::uint32_t>(0u));
  1784. // Extract the data.
  1785. // First get the digits to the left of the decimal point...
  1786. data[0u] = boost::lexical_cast<boost::uint32_t>(str.substr(static_cast<std::size_t>(0u), pos));
  1787. // ...then get the remaining digits to the right of the decimal point.
  1788. const std::string::size_type i_end = ((str.length() - pos_plus_one) / static_cast<std::string::size_type>(cpp_dec_float_elem_digits10));
  1789. for(std::string::size_type i = static_cast<std::string::size_type>(0u); i < i_end; i++)
  1790. {
  1791. const std::string::const_iterator it = str.begin()
  1792. + pos_plus_one
  1793. + (i * static_cast<std::string::size_type>(cpp_dec_float_elem_digits10));
  1794. data[i + 1u] = boost::lexical_cast<boost::uint32_t>(std::string(it, it + static_cast<std::string::size_type>(cpp_dec_float_elem_digits10)));
  1795. }
  1796. // Check for overflow...
  1797. if(exp > cpp_dec_float_max_exp10)
  1798. {
  1799. const bool b_result_is_neg = neg;
  1800. *this = inf();
  1801. if(b_result_is_neg)
  1802. negate();
  1803. }
  1804. // ...and check for underflow.
  1805. if(exp <= cpp_dec_float_min_exp10)
  1806. {
  1807. if(exp == cpp_dec_float_min_exp10)
  1808. {
  1809. // Check for identity with the minimum value.
  1810. cpp_dec_float<Digits10, ExponentType, Allocator> test = *this;
  1811. test.exp = static_cast<ExponentType>(0);
  1812. if(test.isone())
  1813. {
  1814. *this = zero();
  1815. }
  1816. }
  1817. else
  1818. {
  1819. *this = zero();
  1820. }
  1821. }
  1822. }
  1823. catch(const bad_lexical_cast&)
  1824. {
  1825. // Rethrow with better error message:
  1826. std::string msg = "Unable to parse the string \"";
  1827. msg += s;
  1828. msg += "\" as a floating point value.";
  1829. throw std::runtime_error(msg);
  1830. }
  1831. return true;
  1832. }
  1833. template <unsigned Digits10, class ExponentType, class Allocator>
  1834. cpp_dec_float<Digits10, ExponentType, Allocator>::cpp_dec_float(const double mantissa, const ExponentType exponent)
  1835. : data (),
  1836. exp (static_cast<ExponentType>(0)),
  1837. neg (false),
  1838. fpclass (cpp_dec_float_finite),
  1839. prec_elem(cpp_dec_float_elem_number)
  1840. {
  1841. // Create *this cpp_dec_float<Digits10, ExponentType, Allocator> from a given mantissa and exponent.
  1842. // Note: This constructor does not maintain the full precision of double.
  1843. const bool mantissa_is_iszero = (::fabs(mantissa) < ((std::numeric_limits<double>::min)() * (1.0 + std::numeric_limits<double>::epsilon())));
  1844. if(mantissa_is_iszero)
  1845. {
  1846. std::fill(data.begin(), data.end(), static_cast<boost::uint32_t>(0u));
  1847. return;
  1848. }
  1849. const bool b_neg = (mantissa < 0.0);
  1850. double d = ((!b_neg) ? mantissa : -mantissa);
  1851. ExponentType e = exponent;
  1852. while(d > 10.0) { d /= 10.0; ++e; }
  1853. while(d < 1.0) { d *= 10.0; --e; }
  1854. boost::int32_t shift = static_cast<boost::int32_t>(e % static_cast<boost::int32_t>(cpp_dec_float_elem_digits10));
  1855. while(static_cast<boost::int32_t>(shift-- % cpp_dec_float_elem_digits10) != static_cast<boost::int32_t>(0))
  1856. {
  1857. d *= 10.0;
  1858. --e;
  1859. }
  1860. exp = e;
  1861. neg = b_neg;
  1862. std::fill(data.begin(), data.end(), static_cast<boost::uint32_t>(0u));
  1863. static const boost::int32_t digit_ratio = static_cast<boost::int32_t>(static_cast<boost::int32_t>(std::numeric_limits<double>::digits10) / static_cast<boost::int32_t>(cpp_dec_float_elem_digits10));
  1864. static const boost::int32_t digit_loops = static_cast<boost::int32_t>(digit_ratio + static_cast<boost::int32_t>(2));
  1865. for(boost::int32_t i = static_cast<boost::int32_t>(0); i < digit_loops; i++)
  1866. {
  1867. boost::uint32_t n = static_cast<boost::uint32_t>(static_cast<boost::uint64_t>(d));
  1868. data[i] = static_cast<boost::uint32_t>(n);
  1869. d -= static_cast<double>(n);
  1870. d *= static_cast<double>(cpp_dec_float_elem_mask);
  1871. }
  1872. }
  1873. template <unsigned Digits10, class ExponentType, class Allocator>
  1874. cpp_dec_float<Digits10, ExponentType, Allocator>& cpp_dec_float<Digits10, ExponentType, Allocator>::operator= (long double a)
  1875. {
  1876. // Christopher Kormanyos's original code used a cast to long long here, but that fails
  1877. // when long double has more digits than a long long.
  1878. using std::frexp;
  1879. using std::ldexp;
  1880. using std::floor;
  1881. if (a == 0) {
  1882. return *this = zero();
  1883. }
  1884. if (a == 1) {
  1885. return *this = one();
  1886. }
  1887. if((boost::math::isinf)(a))
  1888. {
  1889. return *this = inf();
  1890. }
  1891. if((boost::math::isnan)(a))
  1892. return *this = nan();
  1893. int e;
  1894. long double f, term;
  1895. *this = zero();
  1896. f = frexp(a, &e);
  1897. static const int shift = std::numeric_limits<int>::digits - 1;
  1898. while(f)
  1899. {
  1900. // extract int sized bits from f:
  1901. f = ldexp(f, shift);
  1902. term = floor(f);
  1903. e -= shift;
  1904. *this *= pow2(shift);
  1905. if(term > 0)
  1906. add_unsigned_long_long(static_cast<unsigned>(term));
  1907. else
  1908. sub_unsigned_long_long(static_cast<unsigned>(-term));
  1909. f -= term;
  1910. }
  1911. if(e != 0)
  1912. *this *= pow2(e);
  1913. return *this;
  1914. }
  1915. template <unsigned Digits10, class ExponentType, class Allocator>
  1916. void cpp_dec_float<Digits10, ExponentType, Allocator>::from_unsigned_long_long(const unsigned long long u)
  1917. {
  1918. std::fill(data.begin(), data.end(), static_cast<boost::uint32_t>(0u));
  1919. exp = static_cast<ExponentType>(0);
  1920. neg = false;
  1921. fpclass = cpp_dec_float_finite;
  1922. prec_elem = cpp_dec_float_elem_number;
  1923. std::size_t i =static_cast<std::size_t>(0u);
  1924. unsigned long long uu = u;
  1925. boost::uint32_t temp[(std::numeric_limits<unsigned long long>::digits10 / static_cast<int>(cpp_dec_float_elem_digits10)) + 3] = { static_cast<boost::uint32_t>(0u) };
  1926. while(uu != static_cast<unsigned long long>(0u))
  1927. {
  1928. temp[i] = static_cast<boost::uint32_t>(uu % static_cast<unsigned long long>(cpp_dec_float_elem_mask));
  1929. uu = static_cast<unsigned long long>(uu / static_cast<unsigned long long>(cpp_dec_float_elem_mask));
  1930. ++i;
  1931. }
  1932. if(i > static_cast<std::size_t>(1u))
  1933. {
  1934. exp += static_cast<ExponentType>((i - 1u) * static_cast<std::size_t>(cpp_dec_float_elem_digits10));
  1935. }
  1936. std::reverse(temp, temp + i);
  1937. std::copy(temp, temp + (std::min)(i, static_cast<std::size_t>(cpp_dec_float_elem_number)), data.begin());
  1938. }
  1939. template <unsigned Digits10, class ExponentType, class Allocator>
  1940. boost::uint32_t cpp_dec_float<Digits10, ExponentType, Allocator>::mul_loop_uv(boost::uint32_t* const u, const boost::uint32_t* const v, const boost::int32_t p)
  1941. {
  1942. //
  1943. // There is a limit on how many limbs this algorithm can handle without dropping digits
  1944. // due to overflow in the carry, it is:
  1945. //
  1946. // FLOOR( (2^64 - 1) / (10^8 * 10^8) ) == 1844
  1947. //
  1948. BOOST_STATIC_ASSERT_MSG(cpp_dec_float_elem_number < 1800, "Too many limbs in the data type for the multiplication algorithm - unsupported precision in cpp_dec_float.");
  1949. boost::uint64_t carry = static_cast<boost::uint64_t>(0u);
  1950. for(boost::int32_t j = static_cast<boost::int32_t>(p - 1u); j >= static_cast<boost::int32_t>(0); j--)
  1951. {
  1952. boost::uint64_t sum = carry;
  1953. for(boost::int32_t i = j; i >= static_cast<boost::int32_t>(0); i--)
  1954. {
  1955. sum += static_cast<boost::uint64_t>(u[j - i] * static_cast<boost::uint64_t>(v[i]));
  1956. }
  1957. u[j] = static_cast<boost::uint32_t>(sum % static_cast<boost::uint32_t>(cpp_dec_float_elem_mask));
  1958. carry = static_cast<boost::uint64_t>(sum / static_cast<boost::uint32_t>(cpp_dec_float_elem_mask));
  1959. }
  1960. return static_cast<boost::uint32_t>(carry);
  1961. }
  1962. template <unsigned Digits10, class ExponentType, class Allocator>
  1963. boost::uint32_t cpp_dec_float<Digits10, ExponentType, Allocator>::mul_loop_n(boost::uint32_t* const u, boost::uint32_t n, const boost::int32_t p)
  1964. {
  1965. boost::uint64_t carry = static_cast<boost::uint64_t>(0u);
  1966. // Multiplication loop.
  1967. for(boost::int32_t j = p - 1; j >= static_cast<boost::int32_t>(0); j--)
  1968. {
  1969. const boost::uint64_t t = static_cast<boost::uint64_t>(carry + static_cast<boost::uint64_t>(u[j] * static_cast<boost::uint64_t>(n)));
  1970. carry = static_cast<boost::uint64_t>(t / static_cast<boost::uint32_t>(cpp_dec_float_elem_mask));
  1971. u[j] = static_cast<boost::uint32_t>(t - static_cast<boost::uint64_t>(static_cast<boost::uint32_t>(cpp_dec_float_elem_mask) * static_cast<boost::uint64_t>(carry)));
  1972. }
  1973. return static_cast<boost::uint32_t>(carry);
  1974. }
  1975. template <unsigned Digits10, class ExponentType, class Allocator>
  1976. boost::uint32_t cpp_dec_float<Digits10, ExponentType, Allocator>::div_loop_n(boost::uint32_t* const u, boost::uint32_t n, const boost::int32_t p)
  1977. {
  1978. boost::uint64_t prev = static_cast<boost::uint64_t>(0u);
  1979. for(boost::int32_t j = static_cast<boost::int32_t>(0); j < p; j++)
  1980. {
  1981. const boost::uint64_t t = static_cast<boost::uint64_t>(u[j] + static_cast<boost::uint64_t>(prev * static_cast<boost::uint32_t>(cpp_dec_float_elem_mask)));
  1982. u[j] = static_cast<boost::uint32_t>(t / n);
  1983. prev = static_cast<boost::uint64_t>(t - static_cast<boost::uint64_t>(n * static_cast<boost::uint64_t>(u[j])));
  1984. }
  1985. return static_cast<boost::uint32_t>(prev);
  1986. }
  1987. template <unsigned Digits10, class ExponentType, class Allocator>
  1988. cpp_dec_float<Digits10, ExponentType, Allocator> cpp_dec_float<Digits10, ExponentType, Allocator>::pow2(const long long p)
  1989. {
  1990. // Create a static const table of p^2 for -128 < p < +128.
  1991. // Note: The size of this table must be odd-numbered and
  1992. // symmetric about 0.
  1993. init.do_nothing();
  1994. static const boost::array<cpp_dec_float<Digits10, ExponentType, Allocator>, 255u> p2_data =
  1995. {{
  1996. cpp_dec_float("5.877471754111437539843682686111228389093327783860437607543758531392086297273635864257812500000000000e-39"),
  1997. cpp_dec_float("1.175494350822287507968736537222245677818665556772087521508751706278417259454727172851562500000000000e-38"),
  1998. cpp_dec_float("2.350988701644575015937473074444491355637331113544175043017503412556834518909454345703125000000000000e-38"),
  1999. cpp_dec_float("4.701977403289150031874946148888982711274662227088350086035006825113669037818908691406250000000000000e-38"),
  2000. cpp_dec_float("9.403954806578300063749892297777965422549324454176700172070013650227338075637817382812500000000000000e-38"),
  2001. cpp_dec_float("1.880790961315660012749978459555593084509864890835340034414002730045467615127563476562500000000000000e-37"),
  2002. cpp_dec_float("3.761581922631320025499956919111186169019729781670680068828005460090935230255126953125000000000000000e-37"),
  2003. cpp_dec_float("7.523163845262640050999913838222372338039459563341360137656010920181870460510253906250000000000000000e-37"),
  2004. cpp_dec_float("1.504632769052528010199982767644474467607891912668272027531202184036374092102050781250000000000000000e-36"),
  2005. cpp_dec_float("3.009265538105056020399965535288948935215783825336544055062404368072748184204101562500000000000000000e-36"),
  2006. cpp_dec_float("6.018531076210112040799931070577897870431567650673088110124808736145496368408203125000000000000000000e-36"),
  2007. cpp_dec_float("1.203706215242022408159986214115579574086313530134617622024961747229099273681640625000000000000000000e-35"),
  2008. cpp_dec_float("2.407412430484044816319972428231159148172627060269235244049923494458198547363281250000000000000000000e-35"),
  2009. cpp_dec_float("4.814824860968089632639944856462318296345254120538470488099846988916397094726562500000000000000000000e-35"),
  2010. cpp_dec_float("9.629649721936179265279889712924636592690508241076940976199693977832794189453125000000000000000000000e-35"),
  2011. cpp_dec_float("1.925929944387235853055977942584927318538101648215388195239938795566558837890625000000000000000000000e-34"),
  2012. cpp_dec_float("3.851859888774471706111955885169854637076203296430776390479877591133117675781250000000000000000000000e-34"),
  2013. cpp_dec_float("7.703719777548943412223911770339709274152406592861552780959755182266235351562500000000000000000000000e-34"),
  2014. cpp_dec_float("1.540743955509788682444782354067941854830481318572310556191951036453247070312500000000000000000000000e-33"),
  2015. cpp_dec_float("3.081487911019577364889564708135883709660962637144621112383902072906494140625000000000000000000000000e-33"),
  2016. cpp_dec_float("6.162975822039154729779129416271767419321925274289242224767804145812988281250000000000000000000000000e-33"),
  2017. cpp_dec_float("1.232595164407830945955825883254353483864385054857848444953560829162597656250000000000000000000000000e-32"),
  2018. cpp_dec_float("2.465190328815661891911651766508706967728770109715696889907121658325195312500000000000000000000000000e-32"),
  2019. cpp_dec_float("4.930380657631323783823303533017413935457540219431393779814243316650390625000000000000000000000000000e-32"),
  2020. cpp_dec_float("9.860761315262647567646607066034827870915080438862787559628486633300781250000000000000000000000000000e-32"),
  2021. cpp_dec_float("1.972152263052529513529321413206965574183016087772557511925697326660156250000000000000000000000000000e-31"),
  2022. cpp_dec_float("3.944304526105059027058642826413931148366032175545115023851394653320312500000000000000000000000000000e-31"),
  2023. cpp_dec_float("7.888609052210118054117285652827862296732064351090230047702789306640625000000000000000000000000000000e-31"),
  2024. cpp_dec_float("1.577721810442023610823457130565572459346412870218046009540557861328125000000000000000000000000000000e-30"),
  2025. cpp_dec_float("3.155443620884047221646914261131144918692825740436092019081115722656250000000000000000000000000000000e-30"),
  2026. cpp_dec_float("6.310887241768094443293828522262289837385651480872184038162231445312500000000000000000000000000000000e-30"),
  2027. cpp_dec_float("1.262177448353618888658765704452457967477130296174436807632446289062500000000000000000000000000000000e-29"),
  2028. cpp_dec_float("2.524354896707237777317531408904915934954260592348873615264892578125000000000000000000000000000000000e-29"),
  2029. cpp_dec_float("5.048709793414475554635062817809831869908521184697747230529785156250000000000000000000000000000000000e-29"),
  2030. cpp_dec_float("1.009741958682895110927012563561966373981704236939549446105957031250000000000000000000000000000000000e-28"),
  2031. cpp_dec_float("2.019483917365790221854025127123932747963408473879098892211914062500000000000000000000000000000000000e-28"),
  2032. cpp_dec_float("4.038967834731580443708050254247865495926816947758197784423828125000000000000000000000000000000000000e-28"),
  2033. cpp_dec_float("8.077935669463160887416100508495730991853633895516395568847656250000000000000000000000000000000000000e-28"),
  2034. cpp_dec_float("1.615587133892632177483220101699146198370726779103279113769531250000000000000000000000000000000000000e-27"),
  2035. cpp_dec_float("3.231174267785264354966440203398292396741453558206558227539062500000000000000000000000000000000000000e-27"),
  2036. cpp_dec_float("6.462348535570528709932880406796584793482907116413116455078125000000000000000000000000000000000000000e-27"),
  2037. cpp_dec_float("1.292469707114105741986576081359316958696581423282623291015625000000000000000000000000000000000000000e-26"),
  2038. cpp_dec_float("2.584939414228211483973152162718633917393162846565246582031250000000000000000000000000000000000000000e-26"),
  2039. cpp_dec_float("5.169878828456422967946304325437267834786325693130493164062500000000000000000000000000000000000000000e-26"),
  2040. cpp_dec_float("1.033975765691284593589260865087453566957265138626098632812500000000000000000000000000000000000000000e-25"),
  2041. cpp_dec_float("2.067951531382569187178521730174907133914530277252197265625000000000000000000000000000000000000000000e-25"),
  2042. cpp_dec_float("4.135903062765138374357043460349814267829060554504394531250000000000000000000000000000000000000000000e-25"),
  2043. cpp_dec_float("8.271806125530276748714086920699628535658121109008789062500000000000000000000000000000000000000000000e-25"),
  2044. cpp_dec_float("1.654361225106055349742817384139925707131624221801757812500000000000000000000000000000000000000000000e-24"),
  2045. cpp_dec_float("3.308722450212110699485634768279851414263248443603515625000000000000000000000000000000000000000000000e-24"),
  2046. cpp_dec_float("6.617444900424221398971269536559702828526496887207031250000000000000000000000000000000000000000000000e-24"),
  2047. cpp_dec_float("1.323488980084844279794253907311940565705299377441406250000000000000000000000000000000000000000000000e-23"),
  2048. cpp_dec_float("2.646977960169688559588507814623881131410598754882812500000000000000000000000000000000000000000000000e-23"),
  2049. cpp_dec_float("5.293955920339377119177015629247762262821197509765625000000000000000000000000000000000000000000000000e-23"),
  2050. cpp_dec_float("1.058791184067875423835403125849552452564239501953125000000000000000000000000000000000000000000000000e-22"),
  2051. cpp_dec_float("2.117582368135750847670806251699104905128479003906250000000000000000000000000000000000000000000000000e-22"),
  2052. cpp_dec_float("4.235164736271501695341612503398209810256958007812500000000000000000000000000000000000000000000000000e-22"),
  2053. cpp_dec_float("8.470329472543003390683225006796419620513916015625000000000000000000000000000000000000000000000000000e-22"),
  2054. cpp_dec_float("1.694065894508600678136645001359283924102783203125000000000000000000000000000000000000000000000000000e-21"),
  2055. cpp_dec_float("3.388131789017201356273290002718567848205566406250000000000000000000000000000000000000000000000000000e-21"),
  2056. cpp_dec_float("6.776263578034402712546580005437135696411132812500000000000000000000000000000000000000000000000000000e-21"),
  2057. cpp_dec_float("1.355252715606880542509316001087427139282226562500000000000000000000000000000000000000000000000000000e-20"),
  2058. cpp_dec_float("2.710505431213761085018632002174854278564453125000000000000000000000000000000000000000000000000000000e-20"),
  2059. cpp_dec_float("5.421010862427522170037264004349708557128906250000000000000000000000000000000000000000000000000000000e-20"),
  2060. cpp_dec_float("1.084202172485504434007452800869941711425781250000000000000000000000000000000000000000000000000000000e-19"),
  2061. cpp_dec_float("2.168404344971008868014905601739883422851562500000000000000000000000000000000000000000000000000000000e-19"),
  2062. cpp_dec_float("4.336808689942017736029811203479766845703125000000000000000000000000000000000000000000000000000000000e-19"),
  2063. cpp_dec_float("8.673617379884035472059622406959533691406250000000000000000000000000000000000000000000000000000000000e-19"),
  2064. cpp_dec_float("1.734723475976807094411924481391906738281250000000000000000000000000000000000000000000000000000000000e-18"),
  2065. cpp_dec_float("3.469446951953614188823848962783813476562500000000000000000000000000000000000000000000000000000000000e-18"),
  2066. cpp_dec_float("6.938893903907228377647697925567626953125000000000000000000000000000000000000000000000000000000000000e-18"),
  2067. cpp_dec_float("1.387778780781445675529539585113525390625000000000000000000000000000000000000000000000000000000000000e-17"),
  2068. cpp_dec_float("2.775557561562891351059079170227050781250000000000000000000000000000000000000000000000000000000000000e-17"),
  2069. cpp_dec_float("5.551115123125782702118158340454101562500000000000000000000000000000000000000000000000000000000000000e-17"),
  2070. cpp_dec_float("1.110223024625156540423631668090820312500000000000000000000000000000000000000000000000000000000000000e-16"),
  2071. cpp_dec_float("2.220446049250313080847263336181640625000000000000000000000000000000000000000000000000000000000000000e-16"),
  2072. cpp_dec_float("4.440892098500626161694526672363281250000000000000000000000000000000000000000000000000000000000000000e-16"),
  2073. cpp_dec_float("8.881784197001252323389053344726562500000000000000000000000000000000000000000000000000000000000000000e-16"),
  2074. cpp_dec_float("1.776356839400250464677810668945312500000000000000000000000000000000000000000000000000000000000000000e-15"),
  2075. cpp_dec_float("3.552713678800500929355621337890625000000000000000000000000000000000000000000000000000000000000000000e-15"),
  2076. cpp_dec_float("7.105427357601001858711242675781250000000000000000000000000000000000000000000000000000000000000000000e-15"),
  2077. cpp_dec_float("1.421085471520200371742248535156250000000000000000000000000000000000000000000000000000000000000000000e-14"),
  2078. cpp_dec_float("2.842170943040400743484497070312500000000000000000000000000000000000000000000000000000000000000000000e-14"),
  2079. cpp_dec_float("5.684341886080801486968994140625000000000000000000000000000000000000000000000000000000000000000000000e-14"),
  2080. cpp_dec_float("1.136868377216160297393798828125000000000000000000000000000000000000000000000000000000000000000000000e-13"),
  2081. cpp_dec_float("2.273736754432320594787597656250000000000000000000000000000000000000000000000000000000000000000000000e-13"),
  2082. cpp_dec_float("4.547473508864641189575195312500000000000000000000000000000000000000000000000000000000000000000000000e-13"),
  2083. cpp_dec_float("9.094947017729282379150390625000000000000000000000000000000000000000000000000000000000000000000000000e-13"),
  2084. cpp_dec_float("1.818989403545856475830078125000000000000000000000000000000000000000000000000000000000000000000000000e-12"),
  2085. cpp_dec_float("3.637978807091712951660156250000000000000000000000000000000000000000000000000000000000000000000000000e-12"),
  2086. cpp_dec_float("7.275957614183425903320312500000000000000000000000000000000000000000000000000000000000000000000000000e-12"),
  2087. cpp_dec_float("1.455191522836685180664062500000000000000000000000000000000000000000000000000000000000000000000000000e-11"),
  2088. cpp_dec_float("2.910383045673370361328125000000000000000000000000000000000000000000000000000000000000000000000000000e-11"),
  2089. cpp_dec_float("5.820766091346740722656250000000000000000000000000000000000000000000000000000000000000000000000000000e-11"),
  2090. cpp_dec_float("1.164153218269348144531250000000000000000000000000000000000000000000000000000000000000000000000000000e-10"),
  2091. cpp_dec_float("2.328306436538696289062500000000000000000000000000000000000000000000000000000000000000000000000000000e-10"),
  2092. cpp_dec_float("4.656612873077392578125000000000000000000000000000000000000000000000000000000000000000000000000000000e-10"),
  2093. cpp_dec_float("9.313225746154785156250000000000000000000000000000000000000000000000000000000000000000000000000000000e-10"),
  2094. cpp_dec_float("1.862645149230957031250000000000000000000000000000000000000000000000000000000000000000000000000000000e-9"),
  2095. cpp_dec_float("3.725290298461914062500000000000000000000000000000000000000000000000000000000000000000000000000000000e-9"),
  2096. cpp_dec_float("7.450580596923828125000000000000000000000000000000000000000000000000000000000000000000000000000000000e-9"),
  2097. cpp_dec_float("1.490116119384765625000000000000000000000000000000000000000000000000000000000000000000000000000000000e-8"),
  2098. cpp_dec_float("2.980232238769531250000000000000000000000000000000000000000000000000000000000000000000000000000000000e-8"),
  2099. cpp_dec_float("5.960464477539062500000000000000000000000000000000000000000000000000000000000000000000000000000000000e-8"),
  2100. cpp_dec_float("1.192092895507812500000000000000000000000000000000000000000000000000000000000000000000000000000000000e-7"),
  2101. cpp_dec_float("2.384185791015625000000000000000000000000000000000000000000000000000000000000000000000000000000000000e-7"),
  2102. cpp_dec_float("4.768371582031250000000000000000000000000000000000000000000000000000000000000000000000000000000000000e-7"),
  2103. cpp_dec_float("9.536743164062500000000000000000000000000000000000000000000000000000000000000000000000000000000000000e-7"),
  2104. cpp_dec_float("1.907348632812500000000000000000000000000000000000000000000000000000000000000000000000000000000000000e-6"),
  2105. cpp_dec_float("3.814697265625000000000000000000000000000000000000000000000000000000000000000000000000000000000000000e-6"),
  2106. cpp_dec_float("7.629394531250000000000000000000000000000000000000000000000000000000000000000000000000000000000000000e-6"),
  2107. cpp_dec_float("0.000015258789062500000000000000000000000000000000000000000000000000000000000000000000000000000000000000"),
  2108. cpp_dec_float("0.000030517578125000000000000000000000000000000000000000000000000000000000000000000000000000000000000000"),
  2109. cpp_dec_float("0.000061035156250000000000000000000000000000000000000000000000000000000000000000000000000000000000000000"),
  2110. cpp_dec_float("0.000122070312500000000000000000000000000000000000000000000000000000000000000000000000000000000000000000"),
  2111. cpp_dec_float("0.000244140625000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000"),
  2112. cpp_dec_float("0.000488281250000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000"),
  2113. cpp_dec_float("0.000976562500000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000"),
  2114. cpp_dec_float("0.001953125000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000"),
  2115. cpp_dec_float("0.003906250000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000"),
  2116. cpp_dec_float("0.007812500000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000"),
  2117. cpp_dec_float("0.01562500000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000"),
  2118. cpp_dec_float("0.03125000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000"),
  2119. cpp_dec_float("0.06250000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000"),
  2120. cpp_dec_float("0.125"),
  2121. cpp_dec_float("0.25"),
  2122. cpp_dec_float("0.5"),
  2123. one(),
  2124. two(),
  2125. cpp_dec_float(static_cast<unsigned long long>(4)),
  2126. cpp_dec_float(static_cast<unsigned long long>(8)),
  2127. cpp_dec_float(static_cast<unsigned long long>(16)),
  2128. cpp_dec_float(static_cast<unsigned long long>(32)),
  2129. cpp_dec_float(static_cast<unsigned long long>(64)),
  2130. cpp_dec_float(static_cast<unsigned long long>(128)),
  2131. cpp_dec_float(static_cast<unsigned long long>(256)),
  2132. cpp_dec_float(static_cast<unsigned long long>(512)),
  2133. cpp_dec_float(static_cast<unsigned long long>(1024)),
  2134. cpp_dec_float(static_cast<unsigned long long>(2048)),
  2135. cpp_dec_float(static_cast<unsigned long long>(4096)),
  2136. cpp_dec_float(static_cast<unsigned long long>(8192)),
  2137. cpp_dec_float(static_cast<unsigned long long>(16384)),
  2138. cpp_dec_float(static_cast<unsigned long long>(32768)),
  2139. cpp_dec_float(static_cast<unsigned long long>(65536)),
  2140. cpp_dec_float(static_cast<unsigned long long>(131072)),
  2141. cpp_dec_float(static_cast<unsigned long long>(262144)),
  2142. cpp_dec_float(static_cast<unsigned long long>(524288)),
  2143. cpp_dec_float(static_cast<boost::uint64_t>(1uL << 20u)),
  2144. cpp_dec_float(static_cast<boost::uint64_t>(1uL << 21u)),
  2145. cpp_dec_float(static_cast<boost::uint64_t>(1uL << 22u)),
  2146. cpp_dec_float(static_cast<boost::uint64_t>(1uL << 23u)),
  2147. cpp_dec_float(static_cast<boost::uint64_t>(1uL << 24u)),
  2148. cpp_dec_float(static_cast<boost::uint64_t>(1uL << 25u)),
  2149. cpp_dec_float(static_cast<boost::uint64_t>(1uL << 26u)),
  2150. cpp_dec_float(static_cast<boost::uint64_t>(1uL << 27u)),
  2151. cpp_dec_float(static_cast<boost::uint64_t>(1uL << 28u)),
  2152. cpp_dec_float(static_cast<boost::uint64_t>(1uL << 29u)),
  2153. cpp_dec_float(static_cast<boost::uint64_t>(1uL << 30u)),
  2154. cpp_dec_float(static_cast<boost::uint64_t>(1uL << 31u)),
  2155. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 32u)),
  2156. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 33u)),
  2157. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 34u)),
  2158. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 35u)),
  2159. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 36u)),
  2160. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 37u)),
  2161. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 38u)),
  2162. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 39u)),
  2163. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 40u)),
  2164. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 41u)),
  2165. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 42u)),
  2166. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 43u)),
  2167. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 44u)),
  2168. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 45u)),
  2169. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 46u)),
  2170. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 47u)),
  2171. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 48u)),
  2172. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 49u)),
  2173. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 50u)),
  2174. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 51u)),
  2175. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 52u)),
  2176. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 53u)),
  2177. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 54u)),
  2178. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 55u)),
  2179. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 56u)),
  2180. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 57u)),
  2181. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 58u)),
  2182. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 59u)),
  2183. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 60u)),
  2184. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 61u)),
  2185. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 62u)),
  2186. cpp_dec_float(static_cast<boost::uint64_t>(1uLL << 63u)),
  2187. cpp_dec_float("1.844674407370955161600000000000000000000000000000000000000000000000000000000000000000000000000000000e19"),
  2188. cpp_dec_float("3.689348814741910323200000000000000000000000000000000000000000000000000000000000000000000000000000000e19"),
  2189. cpp_dec_float("7.378697629483820646400000000000000000000000000000000000000000000000000000000000000000000000000000000e19"),
  2190. cpp_dec_float("1.475739525896764129280000000000000000000000000000000000000000000000000000000000000000000000000000000e20"),
  2191. cpp_dec_float("2.951479051793528258560000000000000000000000000000000000000000000000000000000000000000000000000000000e20"),
  2192. cpp_dec_float("5.902958103587056517120000000000000000000000000000000000000000000000000000000000000000000000000000000e20"),
  2193. cpp_dec_float("1.180591620717411303424000000000000000000000000000000000000000000000000000000000000000000000000000000e21"),
  2194. cpp_dec_float("2.361183241434822606848000000000000000000000000000000000000000000000000000000000000000000000000000000e21"),
  2195. cpp_dec_float("4.722366482869645213696000000000000000000000000000000000000000000000000000000000000000000000000000000e21"),
  2196. cpp_dec_float("9.444732965739290427392000000000000000000000000000000000000000000000000000000000000000000000000000000e21"),
  2197. cpp_dec_float("1.888946593147858085478400000000000000000000000000000000000000000000000000000000000000000000000000000e22"),
  2198. cpp_dec_float("3.777893186295716170956800000000000000000000000000000000000000000000000000000000000000000000000000000e22"),
  2199. cpp_dec_float("7.555786372591432341913600000000000000000000000000000000000000000000000000000000000000000000000000000e22"),
  2200. cpp_dec_float("1.511157274518286468382720000000000000000000000000000000000000000000000000000000000000000000000000000e23"),
  2201. cpp_dec_float("3.022314549036572936765440000000000000000000000000000000000000000000000000000000000000000000000000000e23"),
  2202. cpp_dec_float("6.044629098073145873530880000000000000000000000000000000000000000000000000000000000000000000000000000e23"),
  2203. cpp_dec_float("1.208925819614629174706176000000000000000000000000000000000000000000000000000000000000000000000000000e24"),
  2204. cpp_dec_float("2.417851639229258349412352000000000000000000000000000000000000000000000000000000000000000000000000000e24"),
  2205. cpp_dec_float("4.835703278458516698824704000000000000000000000000000000000000000000000000000000000000000000000000000e24"),
  2206. cpp_dec_float("9.671406556917033397649408000000000000000000000000000000000000000000000000000000000000000000000000000e24"),
  2207. cpp_dec_float("1.934281311383406679529881600000000000000000000000000000000000000000000000000000000000000000000000000e25"),
  2208. cpp_dec_float("3.868562622766813359059763200000000000000000000000000000000000000000000000000000000000000000000000000e25"),
  2209. cpp_dec_float("7.737125245533626718119526400000000000000000000000000000000000000000000000000000000000000000000000000e25"),
  2210. cpp_dec_float("1.547425049106725343623905280000000000000000000000000000000000000000000000000000000000000000000000000e26"),
  2211. cpp_dec_float("3.094850098213450687247810560000000000000000000000000000000000000000000000000000000000000000000000000e26"),
  2212. cpp_dec_float("6.189700196426901374495621120000000000000000000000000000000000000000000000000000000000000000000000000e26"),
  2213. cpp_dec_float("1.237940039285380274899124224000000000000000000000000000000000000000000000000000000000000000000000000e27"),
  2214. cpp_dec_float("2.475880078570760549798248448000000000000000000000000000000000000000000000000000000000000000000000000e27"),
  2215. cpp_dec_float("4.951760157141521099596496896000000000000000000000000000000000000000000000000000000000000000000000000e27"),
  2216. cpp_dec_float("9.903520314283042199192993792000000000000000000000000000000000000000000000000000000000000000000000000e27"),
  2217. cpp_dec_float("1.980704062856608439838598758400000000000000000000000000000000000000000000000000000000000000000000000e28"),
  2218. cpp_dec_float("3.961408125713216879677197516800000000000000000000000000000000000000000000000000000000000000000000000e28"),
  2219. cpp_dec_float("7.922816251426433759354395033600000000000000000000000000000000000000000000000000000000000000000000000e28"),
  2220. cpp_dec_float("1.584563250285286751870879006720000000000000000000000000000000000000000000000000000000000000000000000e29"),
  2221. cpp_dec_float("3.169126500570573503741758013440000000000000000000000000000000000000000000000000000000000000000000000e29"),
  2222. cpp_dec_float("6.338253001141147007483516026880000000000000000000000000000000000000000000000000000000000000000000000e29"),
  2223. cpp_dec_float("1.267650600228229401496703205376000000000000000000000000000000000000000000000000000000000000000000000e30"),
  2224. cpp_dec_float("2.535301200456458802993406410752000000000000000000000000000000000000000000000000000000000000000000000e30"),
  2225. cpp_dec_float("5.070602400912917605986812821504000000000000000000000000000000000000000000000000000000000000000000000e30"),
  2226. cpp_dec_float("1.014120480182583521197362564300800000000000000000000000000000000000000000000000000000000000000000000e31"),
  2227. cpp_dec_float("2.028240960365167042394725128601600000000000000000000000000000000000000000000000000000000000000000000e31"),
  2228. cpp_dec_float("4.056481920730334084789450257203200000000000000000000000000000000000000000000000000000000000000000000e31"),
  2229. cpp_dec_float("8.112963841460668169578900514406400000000000000000000000000000000000000000000000000000000000000000000e31"),
  2230. cpp_dec_float("1.622592768292133633915780102881280000000000000000000000000000000000000000000000000000000000000000000e32"),
  2231. cpp_dec_float("3.245185536584267267831560205762560000000000000000000000000000000000000000000000000000000000000000000e32"),
  2232. cpp_dec_float("6.490371073168534535663120411525120000000000000000000000000000000000000000000000000000000000000000000e32"),
  2233. cpp_dec_float("1.298074214633706907132624082305024000000000000000000000000000000000000000000000000000000000000000000e33"),
  2234. cpp_dec_float("2.596148429267413814265248164610048000000000000000000000000000000000000000000000000000000000000000000e33"),
  2235. cpp_dec_float("5.192296858534827628530496329220096000000000000000000000000000000000000000000000000000000000000000000e33"),
  2236. cpp_dec_float("1.038459371706965525706099265844019200000000000000000000000000000000000000000000000000000000000000000e34"),
  2237. cpp_dec_float("2.076918743413931051412198531688038400000000000000000000000000000000000000000000000000000000000000000e34"),
  2238. cpp_dec_float("4.153837486827862102824397063376076800000000000000000000000000000000000000000000000000000000000000000e34"),
  2239. cpp_dec_float("8.307674973655724205648794126752153600000000000000000000000000000000000000000000000000000000000000000e34"),
  2240. cpp_dec_float("1.661534994731144841129758825350430720000000000000000000000000000000000000000000000000000000000000000e35"),
  2241. cpp_dec_float("3.323069989462289682259517650700861440000000000000000000000000000000000000000000000000000000000000000e35"),
  2242. cpp_dec_float("6.646139978924579364519035301401722880000000000000000000000000000000000000000000000000000000000000000e35"),
  2243. cpp_dec_float("1.329227995784915872903807060280344576000000000000000000000000000000000000000000000000000000000000000e36"),
  2244. cpp_dec_float("2.658455991569831745807614120560689152000000000000000000000000000000000000000000000000000000000000000e36"),
  2245. cpp_dec_float("5.316911983139663491615228241121378304000000000000000000000000000000000000000000000000000000000000000e36"),
  2246. cpp_dec_float("1.063382396627932698323045648224275660800000000000000000000000000000000000000000000000000000000000000e37"),
  2247. cpp_dec_float("2.126764793255865396646091296448551321600000000000000000000000000000000000000000000000000000000000000e37"),
  2248. cpp_dec_float("4.253529586511730793292182592897102643200000000000000000000000000000000000000000000000000000000000000e37"),
  2249. cpp_dec_float("8.507059173023461586584365185794205286400000000000000000000000000000000000000000000000000000000000000e37"),
  2250. cpp_dec_float("1.701411834604692317316873037158841057280000000000000000000000000000000000000000000000000000000000000e38")
  2251. }};
  2252. if((p > static_cast<long long>(-128)) && (p < static_cast<long long>(+128)))
  2253. {
  2254. return p2_data[static_cast<std::size_t>(p + ((p2_data.size() - 1u) / 2u))];
  2255. }
  2256. else
  2257. {
  2258. // Compute and return 2^p.
  2259. if(p < static_cast<long long>(0))
  2260. {
  2261. return pow2(static_cast<long long>(-p)).calculate_inv();
  2262. }
  2263. else
  2264. {
  2265. cpp_dec_float<Digits10, ExponentType, Allocator> t;
  2266. default_ops::detail::pow_imp(t, two(), p, mpl::true_());
  2267. return t;
  2268. }
  2269. }
  2270. }
  2271. template <unsigned Digits10, class ExponentType, class Allocator>
  2272. inline void eval_add(cpp_dec_float<Digits10, ExponentType, Allocator>& result, const cpp_dec_float<Digits10, ExponentType, Allocator>& o)
  2273. {
  2274. result += o;
  2275. }
  2276. template <unsigned Digits10, class ExponentType, class Allocator>
  2277. inline void eval_subtract(cpp_dec_float<Digits10, ExponentType, Allocator>& result, const cpp_dec_float<Digits10, ExponentType, Allocator>& o)
  2278. {
  2279. result -= o;
  2280. }
  2281. template <unsigned Digits10, class ExponentType, class Allocator>
  2282. inline void eval_multiply(cpp_dec_float<Digits10, ExponentType, Allocator>& result, const cpp_dec_float<Digits10, ExponentType, Allocator>& o)
  2283. {
  2284. result *= o;
  2285. }
  2286. template <unsigned Digits10, class ExponentType, class Allocator>
  2287. inline void eval_divide(cpp_dec_float<Digits10, ExponentType, Allocator>& result, const cpp_dec_float<Digits10, ExponentType, Allocator>& o)
  2288. {
  2289. result /= o;
  2290. }
  2291. template <unsigned Digits10, class ExponentType, class Allocator>
  2292. inline void eval_add(cpp_dec_float<Digits10, ExponentType, Allocator>& result, const unsigned long long& o)
  2293. {
  2294. result.add_unsigned_long_long(o);
  2295. }
  2296. template <unsigned Digits10, class ExponentType, class Allocator>
  2297. inline void eval_subtract(cpp_dec_float<Digits10, ExponentType, Allocator>& result, const unsigned long long& o)
  2298. {
  2299. result.sub_unsigned_long_long(o);
  2300. }
  2301. template <unsigned Digits10, class ExponentType, class Allocator>
  2302. inline void eval_multiply(cpp_dec_float<Digits10, ExponentType, Allocator>& result, const unsigned long long& o)
  2303. {
  2304. result.mul_unsigned_long_long(o);
  2305. }
  2306. template <unsigned Digits10, class ExponentType, class Allocator>
  2307. inline void eval_divide(cpp_dec_float<Digits10, ExponentType, Allocator>& result, const unsigned long long& o)
  2308. {
  2309. result.div_unsigned_long_long(o);
  2310. }
  2311. template <unsigned Digits10, class ExponentType, class Allocator>
  2312. inline void eval_add(cpp_dec_float<Digits10, ExponentType, Allocator>& result, long long o)
  2313. {
  2314. if(o < 0)
  2315. result.sub_unsigned_long_long(-o);
  2316. else
  2317. result.add_unsigned_long_long(o);
  2318. }
  2319. template <unsigned Digits10, class ExponentType, class Allocator>
  2320. inline void eval_subtract(cpp_dec_float<Digits10, ExponentType, Allocator>& result, long long o)
  2321. {
  2322. if(o < 0)
  2323. result.add_unsigned_long_long(-o);
  2324. else
  2325. result.sub_unsigned_long_long(o);
  2326. }
  2327. template <unsigned Digits10, class ExponentType, class Allocator>
  2328. inline void eval_multiply(cpp_dec_float<Digits10, ExponentType, Allocator>& result, long long o)
  2329. {
  2330. if(o < 0)
  2331. {
  2332. result.mul_unsigned_long_long(-o);
  2333. result.negate();
  2334. }
  2335. else
  2336. result.mul_unsigned_long_long(o);
  2337. }
  2338. template <unsigned Digits10, class ExponentType, class Allocator>
  2339. inline void eval_divide(cpp_dec_float<Digits10, ExponentType, Allocator>& result, long long o)
  2340. {
  2341. if(o < 0)
  2342. {
  2343. result.div_unsigned_long_long(-o);
  2344. result.negate();
  2345. }
  2346. else
  2347. result.div_unsigned_long_long(o);
  2348. }
  2349. template <unsigned Digits10, class ExponentType, class Allocator>
  2350. inline void eval_convert_to(unsigned long long* result, const cpp_dec_float<Digits10, ExponentType, Allocator>& val)
  2351. {
  2352. *result = val.extract_unsigned_long_long();
  2353. }
  2354. template <unsigned Digits10, class ExponentType, class Allocator>
  2355. inline void eval_convert_to(long long* result, const cpp_dec_float<Digits10, ExponentType, Allocator>& val)
  2356. {
  2357. *result = val.extract_signed_long_long();
  2358. }
  2359. template <unsigned Digits10, class ExponentType, class Allocator>
  2360. inline void eval_convert_to(long double* result, cpp_dec_float<Digits10, ExponentType, Allocator>& val)
  2361. {
  2362. *result = val.extract_long_double();
  2363. }
  2364. //
  2365. // Non member function support:
  2366. //
  2367. template <unsigned Digits10, class ExponentType, class Allocator>
  2368. inline int eval_fpclassify(const cpp_dec_float<Digits10, ExponentType, Allocator>& x)
  2369. {
  2370. if((x.isinf)())
  2371. return FP_INFINITE;
  2372. if((x.isnan)())
  2373. return FP_NAN;
  2374. if(x.iszero())
  2375. return FP_ZERO;
  2376. return FP_NORMAL;
  2377. }
  2378. template <unsigned Digits10, class ExponentType, class Allocator>
  2379. inline void eval_abs(cpp_dec_float<Digits10, ExponentType, Allocator>& result, const cpp_dec_float<Digits10, ExponentType, Allocator>& x)
  2380. {
  2381. result = x;
  2382. if(x.isneg())
  2383. result.negate();
  2384. }
  2385. template <unsigned Digits10, class ExponentType, class Allocator>
  2386. inline void eval_fabs(cpp_dec_float<Digits10, ExponentType, Allocator>& result, const cpp_dec_float<Digits10, ExponentType, Allocator>& x)
  2387. {
  2388. result = x;
  2389. if(x.isneg())
  2390. result.negate();
  2391. }
  2392. template <unsigned Digits10, class ExponentType, class Allocator>
  2393. inline void eval_sqrt(cpp_dec_float<Digits10, ExponentType, Allocator>& result, const cpp_dec_float<Digits10, ExponentType, Allocator>& x)
  2394. {
  2395. result = x;
  2396. result.calculate_sqrt();
  2397. }
  2398. template <unsigned Digits10, class ExponentType, class Allocator>
  2399. inline void eval_floor(cpp_dec_float<Digits10, ExponentType, Allocator>& result, const cpp_dec_float<Digits10, ExponentType, Allocator>& x)
  2400. {
  2401. result = x;
  2402. if(!(x.isfinite)() || x.isint())
  2403. {
  2404. return;
  2405. }
  2406. if(x.isneg())
  2407. result -= cpp_dec_float<Digits10, ExponentType, Allocator>::one();
  2408. result = result.extract_integer_part();
  2409. }
  2410. template <unsigned Digits10, class ExponentType, class Allocator>
  2411. inline void eval_ceil(cpp_dec_float<Digits10, ExponentType, Allocator>& result, const cpp_dec_float<Digits10, ExponentType, Allocator>& x)
  2412. {
  2413. result = x;
  2414. if(!(x.isfinite)() || x.isint())
  2415. {
  2416. return;
  2417. }
  2418. if(!x.isneg())
  2419. result += cpp_dec_float<Digits10, ExponentType, Allocator>::one();
  2420. result = result.extract_integer_part();
  2421. }
  2422. template <unsigned Digits10, class ExponentType, class Allocator>
  2423. inline void eval_trunc(cpp_dec_float<Digits10, ExponentType, Allocator>& result, const cpp_dec_float<Digits10, ExponentType, Allocator>& x)
  2424. {
  2425. if(!(x.isfinite)())
  2426. {
  2427. result = boost::math::policies::raise_rounding_error("boost::multiprecision::trunc<%1%>(%1%)", 0, number<cpp_dec_float<Digits10, ExponentType, Allocator> >(x), number<cpp_dec_float<Digits10, ExponentType, Allocator> >(x), boost::math::policies::policy<>()).backend();
  2428. return;
  2429. }
  2430. else if(x.isint())
  2431. {
  2432. result = x;
  2433. return;
  2434. }
  2435. result = x.extract_integer_part();
  2436. }
  2437. template <unsigned Digits10, class ExponentType, class Allocator, class ArgType>
  2438. inline void eval_ldexp(cpp_dec_float<Digits10, ExponentType, Allocator>& result, const cpp_dec_float<Digits10, ExponentType, Allocator>& x, ArgType e)
  2439. {
  2440. const long long the_exp = static_cast<long long>(e);
  2441. if((the_exp > (std::numeric_limits<ExponentType>::max)()) || (the_exp < (std::numeric_limits<ExponentType>::min)()))
  2442. BOOST_THROW_EXCEPTION(std::runtime_error(std::string("Exponent value is out of range.")));
  2443. result = x;
  2444. if ((the_exp > static_cast<long long>(-std::numeric_limits<long long>::digits)) && (the_exp < static_cast<long long>(0)))
  2445. result.div_unsigned_long_long(1ULL << static_cast<long long>(-the_exp));
  2446. else if((the_exp < static_cast<long long>( std::numeric_limits<long long>::digits)) && (the_exp > static_cast<long long>(0)))
  2447. result.mul_unsigned_long_long(1ULL << the_exp);
  2448. else if(the_exp != static_cast<long long>(0))
  2449. result *= cpp_dec_float<Digits10, ExponentType, Allocator>::pow2(e);
  2450. }
  2451. template <unsigned Digits10, class ExponentType, class Allocator>
  2452. inline void eval_frexp(cpp_dec_float<Digits10, ExponentType, Allocator>& result, const cpp_dec_float<Digits10, ExponentType, Allocator>& x, ExponentType* e)
  2453. {
  2454. result = x;
  2455. if(result.isneg())
  2456. result.negate();
  2457. if(result.iszero())
  2458. {
  2459. *e = 0;
  2460. return;
  2461. }
  2462. ExponentType t = result.order();
  2463. BOOST_MP_USING_ABS
  2464. if(abs(t) < ((std::numeric_limits<ExponentType>::max)() / 1000))
  2465. {
  2466. t *= 1000;
  2467. t /= 301;
  2468. }
  2469. else
  2470. {
  2471. t /= 301;
  2472. t *= 1000;
  2473. }
  2474. result *= cpp_dec_float<Digits10, ExponentType, Allocator>::pow2(-t);
  2475. if(result.iszero() || (result.isinf)() || (result.isnan)())
  2476. {
  2477. // pow2 overflowed, slip the calculation up:
  2478. result = x;
  2479. if(result.isneg())
  2480. result.negate();
  2481. t /= 2;
  2482. result *= cpp_dec_float<Digits10, ExponentType, Allocator>::pow2(-t);
  2483. }
  2484. BOOST_MP_USING_ABS
  2485. if(abs(result.order()) > 5)
  2486. {
  2487. // If our first estimate doesn't get close enough then try recursion until we do:
  2488. ExponentType e2;
  2489. cpp_dec_float<Digits10, ExponentType, Allocator> r2;
  2490. eval_frexp(r2, result, &e2);
  2491. // overflow protection:
  2492. if((t > 0) && (e2 > 0) && (t > (std::numeric_limits<ExponentType>::max)() - e2))
  2493. BOOST_THROW_EXCEPTION(std::runtime_error("Exponent is too large to be represented as a power of 2."));
  2494. if((t < 0) && (e2 < 0) && (t < (std::numeric_limits<ExponentType>::min)() - e2))
  2495. BOOST_THROW_EXCEPTION(std::runtime_error("Exponent is too large to be represented as a power of 2."));
  2496. t += e2;
  2497. result = r2;
  2498. }
  2499. while(result.compare(cpp_dec_float<Digits10, ExponentType, Allocator>::one()) >= 0)
  2500. {
  2501. result /= cpp_dec_float<Digits10, ExponentType, Allocator>::two();
  2502. ++t;
  2503. }
  2504. while(result.compare(cpp_dec_float<Digits10, ExponentType, Allocator>::half()) < 0)
  2505. {
  2506. result *= cpp_dec_float<Digits10, ExponentType, Allocator>::two();
  2507. --t;
  2508. }
  2509. *e = t;
  2510. if(x.isneg())
  2511. result.negate();
  2512. }
  2513. template <unsigned Digits10, class ExponentType, class Allocator>
  2514. inline typename disable_if<is_same<ExponentType, int> >::type eval_frexp(cpp_dec_float<Digits10, ExponentType, Allocator>& result, const cpp_dec_float<Digits10, ExponentType, Allocator>& x, int* e)
  2515. {
  2516. ExponentType t;
  2517. eval_frexp(result, x, &t);
  2518. if((t > (std::numeric_limits<int>::max)()) || (t < (std::numeric_limits<int>::min)()))
  2519. BOOST_THROW_EXCEPTION(std::runtime_error("Exponent is outside the range of an int"));
  2520. *e = static_cast<int>(t);
  2521. }
  2522. template <unsigned Digits10, class ExponentType, class Allocator>
  2523. inline bool eval_is_zero(const cpp_dec_float<Digits10, ExponentType, Allocator>& val)
  2524. {
  2525. return val.iszero();
  2526. }
  2527. template <unsigned Digits10, class ExponentType, class Allocator>
  2528. inline int eval_get_sign(const cpp_dec_float<Digits10, ExponentType, Allocator>& val)
  2529. {
  2530. return val.iszero() ? 0 : val.isneg() ? -1 : 1;
  2531. }
  2532. } // namespace backends
  2533. using boost::multiprecision::backends::cpp_dec_float;
  2534. typedef number<cpp_dec_float<50> > cpp_dec_float_50;
  2535. typedef number<cpp_dec_float<100> > cpp_dec_float_100;
  2536. #ifdef BOOST_NO_SFINAE_EXPR
  2537. namespace detail{
  2538. template<unsigned D1, class E1, class A1, unsigned D2, class E2, class A2>
  2539. struct is_explicitly_convertible<cpp_dec_float<D1, E1, A1>, cpp_dec_float<D2, E2, A2> > : public mpl::true_ {};
  2540. }
  2541. #endif
  2542. }}
  2543. namespace std
  2544. {
  2545. template <unsigned Digits10, class ExponentType, class Allocator, boost::multiprecision::expression_template_option ExpressionTemplates>
  2546. class numeric_limits<boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> >
  2547. {
  2548. public:
  2549. BOOST_STATIC_CONSTEXPR bool is_specialized = true;
  2550. BOOST_STATIC_CONSTEXPR bool is_signed = true;
  2551. BOOST_STATIC_CONSTEXPR bool is_integer = false;
  2552. BOOST_STATIC_CONSTEXPR bool is_exact = false;
  2553. BOOST_STATIC_CONSTEXPR bool is_bounded = true;
  2554. BOOST_STATIC_CONSTEXPR bool is_modulo = false;
  2555. BOOST_STATIC_CONSTEXPR bool is_iec559 = false;
  2556. BOOST_STATIC_CONSTEXPR int digits = boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>::cpp_dec_float_digits10;
  2557. BOOST_STATIC_CONSTEXPR int digits10 = boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>::cpp_dec_float_digits10;
  2558. BOOST_STATIC_CONSTEXPR int max_digits10 = boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>::cpp_dec_float_total_digits10;
  2559. BOOST_STATIC_CONSTEXPR ExponentType min_exponent = boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>::cpp_dec_float_min_exp; // Type differs from int.
  2560. BOOST_STATIC_CONSTEXPR ExponentType min_exponent10 = boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>::cpp_dec_float_min_exp10; // Type differs from int.
  2561. BOOST_STATIC_CONSTEXPR ExponentType max_exponent = boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>::cpp_dec_float_max_exp; // Type differs from int.
  2562. BOOST_STATIC_CONSTEXPR ExponentType max_exponent10 = boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>::cpp_dec_float_max_exp10; // Type differs from int.
  2563. BOOST_STATIC_CONSTEXPR int radix = boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>::cpp_dec_float_radix;
  2564. BOOST_STATIC_CONSTEXPR std::float_round_style round_style = std::round_indeterminate;
  2565. BOOST_STATIC_CONSTEXPR bool has_infinity = true;
  2566. BOOST_STATIC_CONSTEXPR bool has_quiet_NaN = true;
  2567. BOOST_STATIC_CONSTEXPR bool has_signaling_NaN = false;
  2568. BOOST_STATIC_CONSTEXPR std::float_denorm_style has_denorm = std::denorm_absent;
  2569. BOOST_STATIC_CONSTEXPR bool has_denorm_loss = false;
  2570. BOOST_STATIC_CONSTEXPR bool traps = false;
  2571. BOOST_STATIC_CONSTEXPR bool tinyness_before = false;
  2572. BOOST_STATIC_CONSTEXPR boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> (min) (void) { return (boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>::min)(); }
  2573. BOOST_STATIC_CONSTEXPR boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> (max) (void) { return (boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>::max)(); }
  2574. BOOST_STATIC_CONSTEXPR boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> lowest (void) { return boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>::zero(); }
  2575. BOOST_STATIC_CONSTEXPR boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> epsilon (void) { return boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>::eps(); }
  2576. BOOST_STATIC_CONSTEXPR boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> round_error (void) { return 0.5L; }
  2577. BOOST_STATIC_CONSTEXPR boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> infinity (void) { return boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>::inf(); }
  2578. BOOST_STATIC_CONSTEXPR boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> quiet_NaN (void) { return boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>::nan(); }
  2579. BOOST_STATIC_CONSTEXPR boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> signaling_NaN(void) { return boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>::zero(); }
  2580. BOOST_STATIC_CONSTEXPR boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> denorm_min (void) { return boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>::zero(); }
  2581. };
  2582. #ifndef BOOST_NO_INCLASS_MEMBER_INITIALIZATION
  2583. template <unsigned Digits10, class ExponentType, class Allocator, boost::multiprecision::expression_template_option ExpressionTemplates>
  2584. BOOST_CONSTEXPR_OR_CONST int numeric_limits<boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> >::digits;
  2585. template <unsigned Digits10, class ExponentType, class Allocator, boost::multiprecision::expression_template_option ExpressionTemplates>
  2586. BOOST_CONSTEXPR_OR_CONST int numeric_limits<boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> >::digits10;
  2587. template <unsigned Digits10, class ExponentType, class Allocator, boost::multiprecision::expression_template_option ExpressionTemplates>
  2588. BOOST_CONSTEXPR_OR_CONST int numeric_limits<boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> >::max_digits10;
  2589. template <unsigned Digits10, class ExponentType, class Allocator, boost::multiprecision::expression_template_option ExpressionTemplates>
  2590. BOOST_CONSTEXPR_OR_CONST bool numeric_limits<boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> >::is_signed;
  2591. template <unsigned Digits10, class ExponentType, class Allocator, boost::multiprecision::expression_template_option ExpressionTemplates>
  2592. BOOST_CONSTEXPR_OR_CONST bool numeric_limits<boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> >::is_integer;
  2593. template <unsigned Digits10, class ExponentType, class Allocator, boost::multiprecision::expression_template_option ExpressionTemplates>
  2594. BOOST_CONSTEXPR_OR_CONST bool numeric_limits<boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> >::is_exact;
  2595. template <unsigned Digits10, class ExponentType, class Allocator, boost::multiprecision::expression_template_option ExpressionTemplates>
  2596. BOOST_CONSTEXPR_OR_CONST int numeric_limits<boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> >::radix;
  2597. template <unsigned Digits10, class ExponentType, class Allocator, boost::multiprecision::expression_template_option ExpressionTemplates>
  2598. BOOST_CONSTEXPR_OR_CONST ExponentType numeric_limits<boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> >::min_exponent;
  2599. template <unsigned Digits10, class ExponentType, class Allocator, boost::multiprecision::expression_template_option ExpressionTemplates>
  2600. BOOST_CONSTEXPR_OR_CONST ExponentType numeric_limits<boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> >::min_exponent10;
  2601. template <unsigned Digits10, class ExponentType, class Allocator, boost::multiprecision::expression_template_option ExpressionTemplates>
  2602. BOOST_CONSTEXPR_OR_CONST ExponentType numeric_limits<boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> >::max_exponent;
  2603. template <unsigned Digits10, class ExponentType, class Allocator, boost::multiprecision::expression_template_option ExpressionTemplates>
  2604. BOOST_CONSTEXPR_OR_CONST ExponentType numeric_limits<boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> >::max_exponent10;
  2605. template <unsigned Digits10, class ExponentType, class Allocator, boost::multiprecision::expression_template_option ExpressionTemplates>
  2606. BOOST_CONSTEXPR_OR_CONST bool numeric_limits<boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> >::has_infinity;
  2607. template <unsigned Digits10, class ExponentType, class Allocator, boost::multiprecision::expression_template_option ExpressionTemplates>
  2608. BOOST_CONSTEXPR_OR_CONST bool numeric_limits<boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> >::has_quiet_NaN;
  2609. template <unsigned Digits10, class ExponentType, class Allocator, boost::multiprecision::expression_template_option ExpressionTemplates>
  2610. BOOST_CONSTEXPR_OR_CONST bool numeric_limits<boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> >::has_signaling_NaN;
  2611. template <unsigned Digits10, class ExponentType, class Allocator, boost::multiprecision::expression_template_option ExpressionTemplates>
  2612. BOOST_CONSTEXPR_OR_CONST float_denorm_style numeric_limits<boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> >::has_denorm;
  2613. template <unsigned Digits10, class ExponentType, class Allocator, boost::multiprecision::expression_template_option ExpressionTemplates>
  2614. BOOST_CONSTEXPR_OR_CONST bool numeric_limits<boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> >::has_denorm_loss;
  2615. template <unsigned Digits10, class ExponentType, class Allocator, boost::multiprecision::expression_template_option ExpressionTemplates>
  2616. BOOST_CONSTEXPR_OR_CONST bool numeric_limits<boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> >::is_iec559;
  2617. template <unsigned Digits10, class ExponentType, class Allocator, boost::multiprecision::expression_template_option ExpressionTemplates>
  2618. BOOST_CONSTEXPR_OR_CONST bool numeric_limits<boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> >::is_bounded;
  2619. template <unsigned Digits10, class ExponentType, class Allocator, boost::multiprecision::expression_template_option ExpressionTemplates>
  2620. BOOST_CONSTEXPR_OR_CONST bool numeric_limits<boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> >::is_modulo;
  2621. template <unsigned Digits10, class ExponentType, class Allocator, boost::multiprecision::expression_template_option ExpressionTemplates>
  2622. BOOST_CONSTEXPR_OR_CONST bool numeric_limits<boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> >::traps;
  2623. template <unsigned Digits10, class ExponentType, class Allocator, boost::multiprecision::expression_template_option ExpressionTemplates>
  2624. BOOST_CONSTEXPR_OR_CONST bool numeric_limits<boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> >::tinyness_before;
  2625. template <unsigned Digits10, class ExponentType, class Allocator, boost::multiprecision::expression_template_option ExpressionTemplates>
  2626. BOOST_CONSTEXPR_OR_CONST float_round_style numeric_limits<boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates> >::round_style;
  2627. #endif
  2628. }
  2629. namespace boost{ namespace math{
  2630. namespace policies{
  2631. template <unsigned Digits10, class ExponentType, class Allocator, class Policy, boost::multiprecision::expression_template_option ExpressionTemplates>
  2632. struct precision< boost::multiprecision::number<boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>, ExpressionTemplates>, Policy>
  2633. {
  2634. // Define a local copy of cpp_dec_float_digits10 because it might differ
  2635. // from the template parameter Digits10 for small or large digit counts.
  2636. static const boost::int32_t cpp_dec_float_digits10 = boost::multiprecision::cpp_dec_float<Digits10, ExponentType, Allocator>::cpp_dec_float_digits10;
  2637. typedef typename Policy::precision_type precision_type;
  2638. typedef digits2<((cpp_dec_float_digits10 + 1LL) * 1000LL) / 301LL> digits_2;
  2639. typedef typename mpl::if_c<
  2640. ((digits_2::value <= precision_type::value)
  2641. || (Policy::precision_type::value <= 0)),
  2642. // Default case, full precision for RealType:
  2643. digits_2,
  2644. // User customised precision:
  2645. precision_type
  2646. >::type type;
  2647. };
  2648. } // namespace policies
  2649. }} // namespaces boost::math
  2650. #endif