1234567891011121314151617181920212223242526272829303132333435363738394041424344454647484950515253545556575859606162636465666768697071727374757677787980818283848586878889909192939495969798991001011021031041051061071081091101111121131141151161171181191201211221231241251261271281291301311321331341351361371381391401411421431441451461471481491501511521531541551561571581591601611621631641651661671681691701711721731741751761771781791801811821831841851861871881891901911921931941951961971981992002012022032042052062072082092102112122132142152162172182192202212222232242252262272282292302312322332342352362372382392402412422432442452462472482492502512522532542552562572582592602612622632642652662672682692702712722732742752762772782792802812822832842852862872882892902912922932942952962972982993003013023033043053063073083093103113123133143153163173183193203213223233243253263273283293303313323333343353363373383393403413423433443453463473483493503513523533543553563573583593603613623633643653663673683693703713723733743753763773783793803813823833843853863873883893903913923933943953963973983994004014024034044054064074084094104114124134144154164174184194204214224234244254264274284294304314324334344354364374384394404414424434444454464474484494504514524534544554564574584594604614624634644654664674684694704714724734744754764774784794804814824834844854864874884894904914924934944954964974984995005015025035045055065075085095105115125135145155165175185195205215225235245255265275285295305315325335345355365375385395405415425435445455465475485495505515525535545555565575585595605615625635645655665675685695705715725735745755765775785795805815825835845855865875885895905915925935945955965975985996006016026036046056066076086096106116126136146156166176186196206216226236246256266276286296306316326336346356366376386396406416426436446456466476486496506516526536546556566576586596606616626636646656666676686696706716726736746756766776786796806816826836846856866876886896906916926936946956966976986997007017027037047057067077087097107117127137147157167177187197207217227237247257267277287297307317327337347357367377387397407417427437447457467477487497507517527537547557567577587597607617627637647657667677687697707717727737747757767777787797807817827837847857867877887897907917927937947957967977987998008018028038048058068078088098108118128138148158168178188198208218228238248258268278288298308318328338348358368378388398408418428438448458468478488498508518528538548558568578588598608618628638648658668678688698708718728738748758768778788798808818828838848858868878888898908918928938948958968978988999009019029039049059069079089099109119129139149159169179189199209219229239249259269279289299309319329339349359369379389399409419429439449459469479489499509519529539549559569579589599609619629639649659669679689699709719729739749759769779789799809819829839849859869879889899909919929939949959969979989991000100110021003100410051006100710081009101010111012101310141015101610171018101910201021102210231024102510261027102810291030103110321033103410351036103710381039104010411042104310441045104610471048104910501051105210531054105510561057105810591060106110621063106410651066106710681069107010711072107310741075107610771078107910801081108210831084108510861087108810891090109110921093109410951096109710981099110011011102110311041105110611071108110911101111111211131114111511161117111811191120112111221123112411251126112711281129113011311132113311341135113611371138113911401141114211431144114511461147114811491150115111521153115411551156115711581159116011611162116311641165116611671168116911701171117211731174117511761177117811791180118111821183118411851186118711881189119011911192119311941195119611971198119912001201120212031204120512061207120812091210121112121213121412151216121712181219122012211222122312241225122612271228122912301231123212331234123512361237123812391240124112421243124412451246124712481249125012511252125312541255125612571258125912601261126212631264126512661267126812691270127112721273127412751276127712781279128012811282128312841285128612871288128912901291129212931294129512961297129812991300130113021303130413051306130713081309131013111312131313141315131613171318131913201321132213231324132513261327132813291330133113321333133413351336133713381339134013411342134313441345134613471348134913501351135213531354135513561357135813591360136113621363136413651366136713681369137013711372137313741375137613771378137913801381138213831384138513861387138813891390139113921393139413951396139713981399140014011402140314041405140614071408140914101411141214131414141514161417141814191420142114221423142414251426142714281429143014311432143314341435143614371438143914401441144214431444144514461447144814491450145114521453145414551456145714581459146014611462146314641465146614671468146914701471147214731474147514761477147814791480148114821483148414851486148714881489149014911492149314941495149614971498149915001501150215031504150515061507150815091510151115121513151415151516151715181519152015211522152315241525152615271528152915301531153215331534153515361537153815391540154115421543154415451546154715481549155015511552155315541555155615571558155915601561156215631564156515661567156815691570157115721573157415751576157715781579158015811582158315841585158615871588158915901591159215931594159515961597159815991600160116021603160416051606160716081609161016111612161316141615161616171618161916201621162216231624162516261627162816291630163116321633163416351636163716381639164016411642164316441645164616471648164916501651165216531654165516561657165816591660166116621663166416651666166716681669167016711672167316741675167616771678167916801681168216831684168516861687168816891690169116921693169416951696169716981699170017011702170317041705170617071708170917101711171217131714171517161717171817191720172117221723172417251726172717281729173017311732173317341735173617371738173917401741174217431744174517461747174817491750175117521753175417551756175717581759176017611762176317641765176617671768176917701771177217731774177517761777177817791780178117821783178417851786178717881789179017911792179317941795179617971798179918001801180218031804180518061807180818091810181118121813181418151816181718181819182018211822182318241825182618271828182918301831183218331834183518361837183818391840184118421843184418451846184718481849185018511852185318541855185618571858185918601861186218631864186518661867186818691870187118721873187418751876187718781879188018811882188318841885188618871888188918901891189218931894189518961897189818991900190119021903190419051906190719081909191019111912191319141915191619171918191919201921192219231924192519261927192819291930193119321933193419351936193719381939194019411942194319441945194619471948194919501951195219531954195519561957195819591960196119621963196419651966196719681969197019711972197319741975197619771978197919801981198219831984198519861987198819891990199119921993199419951996199719981999200020012002200320042005200620072008200920102011201220132014201520162017201820192020202120222023202420252026202720282029203020312032203320342035203620372038203920402041204220432044204520462047204820492050205120522053205420552056205720582059206020612062206320642065206620672068206920702071207220732074207520762077207820792080208120822083208420852086208720882089209020912092209320942095209620972098209921002101210221032104210521062107210821092110211121122113211421152116211721182119212021212122212321242125212621272128212921302131213221332134213521362137213821392140214121422143214421452146214721482149215021512152215321542155215621572158215921602161216221632164216521662167216821692170217121722173217421752176217721782179218021812182218321842185218621872188218921902191219221932194219521962197219821992200220122022203220422052206220722082209221022112212221322142215221622172218221922202221222222232224222522262227222822292230223122322233223422352236223722382239224022412242224322442245224622472248224922502251225222532254225522562257225822592260226122622263226422652266226722682269227022712272227322742275227622772278227922802281228222832284228522862287228822892290229122922293229422952296229722982299230023012302230323042305230623072308230923102311231223132314231523162317231823192320232123222323232423252326232723282329233023312332233323342335233623372338233923402341234223432344234523462347234823492350235123522353235423552356235723582359236023612362236323642365236623672368236923702371237223732374237523762377237823792380238123822383238423852386238723882389239023912392239323942395239623972398239924002401240224032404240524062407240824092410241124122413241424152416241724182419242024212422242324242425242624272428242924302431243224332434243524362437243824392440244124422443244424452446244724482449245024512452245324542455245624572458245924602461246224632464246524662467246824692470247124722473247424752476247724782479248024812482248324842485248624872488248924902491249224932494249524962497249824992500250125022503250425052506250725082509251025112512251325142515251625172518251925202521252225232524252525262527252825292530253125322533253425352536253725382539254025412542254325442545254625472548254925502551255225532554255525562557255825592560256125622563256425652566256725682569257025712572257325742575257625772578257925802581258225832584258525862587 |
- /* Copyright (C) 2007-2022 Free Software Foundation, Inc.
- This file is part of GCC.
- GCC is free software; you can redistribute it and/or modify it under
- the terms of the GNU General Public License as published by the Free
- Software Foundation; either version 3, or (at your option) any later
- version.
- GCC is distributed in the hope that it will be useful, but WITHOUT ANY
- WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
- Under Section 7 of GPL version 3, you are granted additional
- permissions described in the GCC Runtime Library Exception, version
- 3.1, as published by the Free Software Foundation.
- You should have received a copy of the GNU General Public License and
- a copy of the GCC Runtime Library Exception along with this program;
- see the files COPYING3 and COPYING.RUNTIME respectively. If not, see
- <http://www.gnu.org/licenses/>. */
- #include "bid_internal.h"
- /*****************************************************************************
- * BID64_to_int32_rnint
- ****************************************************************************/
- #if DECIMAL_CALL_BY_REFERENCE
- void
- bid64_to_int32_rnint (int *pres,
- UINT64 *
- px _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
- _EXC_INFO_PARAM) {
- UINT64 x = *px;
- #else
- int
- bid64_to_int32_rnint (UINT64 x _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
- _EXC_INFO_PARAM) {
- #endif
- int res;
- UINT64 x_sign;
- UINT64 x_exp;
- int exp; // unbiased exponent
- // Note: C1 represents x_significand (UINT64)
- UINT64 tmp64;
- BID_UI64DOUBLE tmp1;
- unsigned int x_nr_bits;
- int q, ind, shift;
- UINT64 C1;
- UINT64 Cstar; // C* represents up to 16 decimal digits ~ 54 bits
- UINT128 fstar;
- UINT128 P128;
- // check for NaN or Infinity
- if ((x & MASK_NAN) == MASK_NAN || (x & MASK_INF) == MASK_INF) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // unpack x
- x_sign = x & MASK_SIGN; // 0 for positive, MASK_SIGN for negative
- // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] =>
- if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) {
- x_exp = (x & MASK_BINARY_EXPONENT2) >> 51; // biased
- C1 = (x & MASK_BINARY_SIG2) | MASK_BINARY_OR2;
- if (C1 > 9999999999999999ull) { // non-canonical
- x_exp = 0;
- C1 = 0;
- }
- } else {
- x_exp = (x & MASK_BINARY_EXPONENT1) >> 53; // biased
- C1 = x & MASK_BINARY_SIG1;
- }
- // check for zeros (possibly from non-canonical values)
- if (C1 == 0x0ull) {
- // x is 0
- res = 0x00000000;
- BID_RETURN (res);
- }
- // x is not special and is not zero
- // q = nr. of decimal digits in x (1 <= q <= 54)
- // determine first the nr. of bits in x
- if (C1 >= 0x0020000000000000ull) { // x >= 2^53
- // split the 64-bit value in two 32-bit halves to avoid rounding errors
- if (C1 >= 0x0000000100000000ull) { // x >= 2^32
- tmp1.d = (double) (C1 >> 32); // exact conversion
- x_nr_bits =
- 33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- } else { // x < 2^32
- tmp1.d = (double) C1; // exact conversion
- x_nr_bits =
- 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- }
- } else { // if x < 2^53
- tmp1.d = (double) C1; // exact conversion
- x_nr_bits =
- 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- }
- q = nr_digits[x_nr_bits - 1].digits;
- if (q == 0) {
- q = nr_digits[x_nr_bits - 1].digits1;
- if (C1 >= nr_digits[x_nr_bits - 1].threshold_lo)
- q++;
- }
- exp = x_exp - 398; // unbiased exponent
- if ((q + exp) > 10) { // x >= 10^10 ~= 2^33.2... (cannot fit in 32 bits)
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- } else if ((q + exp) == 10) { // x = c(0)c(1)...c(9).c(10)...c(q-1)
- // in this case 2^29.89... ~= 10^9 <= x < 10^10 ~= 2^33.2...
- // so x rounded to an integer may or may not fit in a signed 32-bit int
- // the cases that do not fit are identified here; the ones that fit
- // fall through and will be handled with other cases further,
- // under '1 <= q + exp <= 10'
- if (x_sign) { // if n < 0 and q + exp = 10
- // if n < -2^31 - 1/2 then n is too large
- // too large if c(0)c(1)...c(9).c(10)...c(q-1) > 2^31+1/2
- // <=> 0.c(0)c(1)...c(q-1) * 10^11 > 0x500000005, 1<=q<=16
- // <=> C * 10^(11-q) > 0x500000005, 1<=q<=16
- if (q <= 11) {
- // Note: C * 10^(11-q) has 10 or 11 digits; 0x500000005 has 11 digits
- tmp64 = C1 * ten2k64[11 - q]; // C scaled up to 11-digit int
- // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
- if (tmp64 > 0x500000005ull) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- } else { // if (q > 11), i.e. 12 <= q <= 16 and so -15 <= exp <= -2
- // C * 10^(11-q) > 0x500000005 <=>
- // C > 0x500000005 * 10^(q-11) where 1 <= q - 11 <= 5
- // (scale 2^31+1/2 up)
- // Note: 0x500000005*10^(q-11) has q-1 or q digits, where q <= 16
- tmp64 = 0x500000005ull * ten2k64[q - 11];
- if (C1 > tmp64) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- }
- } else { // if n > 0 and q + exp = 10
- // if n >= 2^31 - 1/2 then n is too large
- // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31-1/2
- // <=> 0.c(0)c(1)...c(q-1) * 10^11 >= 0x4fffffffb, 1<=q<=16
- // <=> C * 10^(11-q) >= 0x4fffffffb, 1<=q<=16
- if (q <= 11) {
- // Note: C * 10^(11-q) has 10 or 11 digits; 0x500000005 has 11 digits
- tmp64 = C1 * ten2k64[11 - q]; // C scaled up to 11-digit int
- // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
- if (tmp64 >= 0x4fffffffbull) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- } else { // if (q > 11), i.e. 12 <= q <= 16 and so -15 <= exp <= -2
- // C * 10^(11-q) >= 0x4fffffffb <=>
- // C >= 0x4fffffffb * 10^(q-11) where 1 <= q - 11 <= 5
- // (scale 2^31-1/2 up)
- // Note: 0x4fffffffb*10^(q-11) has q-1 or q digits, where q <= 16
- tmp64 = 0x4fffffffbull * ten2k64[q - 11];
- if (C1 >= tmp64) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- }
- }
- }
- // n is not too large to be converted to int32: -2^31 - 1/2 <= n < 2^31 - 1/2
- // Note: some of the cases tested for above fall through to this point
- if ((q + exp) < 0) { // n = +/-0.0...c(0)c(1)...c(q-1)
- // return 0
- res = 0x00000000;
- BID_RETURN (res);
- } else if ((q + exp) == 0) { // n = +/-0.c(0)c(1)...c(q-1)
- // if 0.c(0)c(1)...c(q-1) <= 0.5 <=> c(0)c(1)...c(q-1) <= 5 * 10^(q-1)
- // res = 0
- // else
- // res = +/-1
- ind = q - 1;
- if (C1 <= midpoint64[ind]) {
- res = 0x00000000; // return 0
- } else if (x_sign) { // n < 0
- res = 0xffffffff; // return -1
- } else { // n > 0
- res = 0x00000001; // return +1
- }
- } else { // if (1 <= q + exp <= 10, 1 <= q <= 16, -15 <= exp <= 9)
- // -2^31-1/2 <= x <= -1 or 1 <= x < 2^31-1/2 so x can be rounded
- // to nearest to a 32-bit signed integer
- if (exp < 0) { // 2 <= q <= 16, -15 <= exp <= -1, 1 <= q + exp <= 10
- ind = -exp; // 1 <= ind <= 15; ind is a synonym for 'x'
- // chop off ind digits from the lower part of C1
- // C1 = C1 + 1/2 * 10^ind where the result C1 fits in 64 bits
- C1 = C1 + midpoint64[ind - 1];
- // calculate C* and f*
- // C* is actually floor(C*) in this case
- // C* and f* need shifting and masking, as shown by
- // shiftright128[] and maskhigh128[]
- // 1 <= x <= 15
- // kx = 10^(-x) = ten2mk64[ind - 1]
- // C* = (C1 + 1/2 * 10^x) * 10^(-x)
- // the approximation of 10^(-x) was rounded up to 54 bits
- __mul_64x64_to_128MACH (P128, C1, ten2mk64[ind - 1]);
- Cstar = P128.w[1];
- fstar.w[1] = P128.w[1] & maskhigh128[ind - 1];
- fstar.w[0] = P128.w[0];
- // the top Ex bits of 10^(-x) are T* = ten2mk128trunc[ind].w[0], e.g.
- // if x=1, T*=ten2mk128trunc[0].w[0]=0x1999999999999999
- // if (0 < f* < 10^(-x)) then the result is a midpoint
- // if floor(C*) is even then C* = floor(C*) - logical right
- // shift; C* has p decimal digits, correct by Prop. 1)
- // else if floor(C*) is odd C* = floor(C*)-1 (logical right
- // shift; C* has p decimal digits, correct by Pr. 1)
- // else
- // C* = floor(C*) (logical right shift; C has p decimal digits,
- // correct by Property 1)
- // n = C* * 10^(e+x)
- // shift right C* by Ex-64 = shiftright128[ind]
- shift = shiftright128[ind - 1]; // 0 <= shift <= 39
- Cstar = Cstar >> shift;
- // if the result was a midpoint it was rounded away from zero, so
- // it will need a correction
- // check for midpoints
- if ((fstar.w[1] == 0) && fstar.w[0]
- && (fstar.w[0] <= ten2mk128trunc[ind - 1].w[1])) {
- // ten2mk128trunc[ind -1].w[1] is identical to
- // ten2mk128[ind -1].w[1]
- // the result is a midpoint; round to nearest
- if (Cstar & 0x01) { // Cstar is odd; MP in [EVEN, ODD]
- // if floor(C*) is odd C = floor(C*) - 1; the result >= 1
- Cstar--; // Cstar is now even
- } // else MP in [ODD, EVEN]
- }
- if (x_sign)
- res = -Cstar;
- else
- res = Cstar;
- } else if (exp == 0) {
- // 1 <= q <= 10
- // res = +/-C (exact)
- if (x_sign)
- res = -C1;
- else
- res = C1;
- } else { // if (exp > 0) => 1 <= exp <= 9, 1 <= q < 9, 2 <= q + exp <= 10
- // res = +/-C * 10^exp (exact)
- if (x_sign)
- res = -C1 * ten2k64[exp];
- else
- res = C1 * ten2k64[exp];
- }
- }
- BID_RETURN (res);
- }
- /*****************************************************************************
- * BID64_to_int32_xrnint
- ****************************************************************************/
- #if DECIMAL_CALL_BY_REFERENCE
- void
- bid64_to_int32_xrnint (int *pres,
- UINT64 *
- px _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
- _EXC_INFO_PARAM) {
- UINT64 x = *px;
- #else
- int
- bid64_to_int32_xrnint (UINT64 x _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
- _EXC_INFO_PARAM) {
- #endif
- int res;
- UINT64 x_sign;
- UINT64 x_exp;
- int exp; // unbiased exponent
- // Note: C1 represents x_significand (UINT64)
- UINT64 tmp64;
- BID_UI64DOUBLE tmp1;
- unsigned int x_nr_bits;
- int q, ind, shift;
- UINT64 C1;
- UINT64 Cstar; // C* represents up to 16 decimal digits ~ 54 bits
- UINT128 fstar;
- UINT128 P128;
- // check for NaN or Infinity
- if ((x & MASK_NAN) == MASK_NAN || (x & MASK_INF) == MASK_INF) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // unpack x
- x_sign = x & MASK_SIGN; // 0 for positive, MASK_SIGN for negative
- // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] =>
- if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) {
- x_exp = (x & MASK_BINARY_EXPONENT2) >> 51; // biased
- C1 = (x & MASK_BINARY_SIG2) | MASK_BINARY_OR2;
- if (C1 > 9999999999999999ull) { // non-canonical
- x_exp = 0;
- C1 = 0;
- }
- } else {
- x_exp = (x & MASK_BINARY_EXPONENT1) >> 53; // biased
- C1 = x & MASK_BINARY_SIG1;
- }
- // check for zeros (possibly from non-canonical values)
- if (C1 == 0x0ull) {
- // x is 0
- res = 0x00000000;
- BID_RETURN (res);
- }
- // x is not special and is not zero
- // q = nr. of decimal digits in x (1 <= q <= 54)
- // determine first the nr. of bits in x
- if (C1 >= 0x0020000000000000ull) { // x >= 2^53
- // split the 64-bit value in two 32-bit halves to avoid rounding errors
- if (C1 >= 0x0000000100000000ull) { // x >= 2^32
- tmp1.d = (double) (C1 >> 32); // exact conversion
- x_nr_bits =
- 33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- } else { // x < 2^32
- tmp1.d = (double) C1; // exact conversion
- x_nr_bits =
- 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- }
- } else { // if x < 2^53
- tmp1.d = (double) C1; // exact conversion
- x_nr_bits =
- 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- }
- q = nr_digits[x_nr_bits - 1].digits;
- if (q == 0) {
- q = nr_digits[x_nr_bits - 1].digits1;
- if (C1 >= nr_digits[x_nr_bits - 1].threshold_lo)
- q++;
- }
- exp = x_exp - 398; // unbiased exponent
- if ((q + exp) > 10) { // x >= 10^10 ~= 2^33.2... (cannot fit in 32 bits)
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- } else if ((q + exp) == 10) { // x = c(0)c(1)...c(9).c(10)...c(q-1)
- // in this case 2^29.89... ~= 10^9 <= x < 10^10 ~= 2^33.2...
- // so x rounded to an integer may or may not fit in a signed 32-bit int
- // the cases that do not fit are identified here; the ones that fit
- // fall through and will be handled with other cases further,
- // under '1 <= q + exp <= 10'
- if (x_sign) { // if n < 0 and q + exp = 10
- // if n < -2^31 - 1/2 then n is too large
- // too large if c(0)c(1)...c(9).c(10)...c(q-1) > 2^31+1/2
- // <=> 0.c(0)c(1)...c(q-1) * 10^11 > 0x500000005, 1<=q<=16
- // <=> C * 10^(11-q) > 0x500000005, 1<=q<=16
- if (q <= 11) {
- // Note: C * 10^(11-q) has 10 or 11 digits; 0x500000005 has 11 digits
- tmp64 = C1 * ten2k64[11 - q]; // C scaled up to 11-digit int
- // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
- if (tmp64 > 0x500000005ull) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- } else { // if (q > 11), i.e. 12 <= q <= 16 and so -15 <= exp <= -2
- // C * 10^(11-q) > 0x500000005 <=>
- // C > 0x500000005 * 10^(q-11) where 1 <= q - 11 <= 5
- // (scale 2^31+1/2 up)
- // Note: 0x500000005*10^(q-11) has q-1 or q digits, where q <= 16
- tmp64 = 0x500000005ull * ten2k64[q - 11];
- if (C1 > tmp64) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- }
- } else { // if n > 0 and q + exp = 10
- // if n >= 2^31 - 1/2 then n is too large
- // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31-1/2
- // <=> 0.c(0)c(1)...c(q-1) * 10^11 >= 0x4fffffffb, 1<=q<=16
- // <=> C * 10^(11-q) >= 0x4fffffffb, 1<=q<=16
- if (q <= 11) {
- // Note: C * 10^(11-q) has 10 or 11 digits; 0x500000005 has 11 digits
- tmp64 = C1 * ten2k64[11 - q]; // C scaled up to 11-digit int
- // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
- if (tmp64 >= 0x4fffffffbull) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- } else { // if (q > 11), i.e. 12 <= q <= 16 and so -15 <= exp <= -2
- // C * 10^(11-q) >= 0x4fffffffb <=>
- // C >= 0x4fffffffb * 10^(q-11) where 1 <= q - 11 <= 5
- // (scale 2^31-1/2 up)
- // Note: 0x4fffffffb*10^(q-11) has q-1 or q digits, where q <= 16
- tmp64 = 0x4fffffffbull * ten2k64[q - 11];
- if (C1 >= tmp64) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- }
- }
- }
- // n is not too large to be converted to int32: -2^31 - 1/2 < n < 2^31 - 1/2
- // Note: some of the cases tested for above fall through to this point
- if ((q + exp) < 0) { // n = +/-0.0...c(0)c(1)...c(q-1)
- // set inexact flag
- *pfpsf |= INEXACT_EXCEPTION;
- // return 0
- res = 0x00000000;
- BID_RETURN (res);
- } else if ((q + exp) == 0) { // n = +/-0.c(0)c(1)...c(q-1)
- // if 0.c(0)c(1)...c(q-1) <= 0.5 <=> c(0)c(1)...c(q-1) <= 5 * 10^(q-1)
- // res = 0
- // else
- // res = +/-1
- ind = q - 1;
- if (C1 <= midpoint64[ind]) {
- res = 0x00000000; // return 0
- } else if (x_sign) { // n < 0
- res = 0xffffffff; // return -1
- } else { // n > 0
- res = 0x00000001; // return +1
- }
- // set inexact flag
- *pfpsf |= INEXACT_EXCEPTION;
- } else { // if (1 <= q + exp <= 10, 1 <= q <= 16, -15 <= exp <= 9)
- // -2^31-1/2 <= x <= -1 or 1 <= x < 2^31-1/2 so x can be rounded
- // to nearest to a 32-bit signed integer
- if (exp < 0) { // 2 <= q <= 16, -15 <= exp <= -1, 1 <= q + exp <= 10
- ind = -exp; // 1 <= ind <= 15; ind is a synonym for 'x'
- // chop off ind digits from the lower part of C1
- // C1 = C1 + 1/2 * 10^ind where the result C1 fits in 64 bits
- C1 = C1 + midpoint64[ind - 1];
- // calculate C* and f*
- // C* is actually floor(C*) in this case
- // C* and f* need shifting and masking, as shown by
- // shiftright128[] and maskhigh128[]
- // 1 <= x <= 15
- // kx = 10^(-x) = ten2mk64[ind - 1]
- // C* = (C1 + 1/2 * 10^x) * 10^(-x)
- // the approximation of 10^(-x) was rounded up to 54 bits
- __mul_64x64_to_128MACH (P128, C1, ten2mk64[ind - 1]);
- Cstar = P128.w[1];
- fstar.w[1] = P128.w[1] & maskhigh128[ind - 1];
- fstar.w[0] = P128.w[0];
- // the top Ex bits of 10^(-x) are T* = ten2mk128trunc[ind].w[0], e.g.
- // if x=1, T*=ten2mk128trunc[0].w[0]=0x1999999999999999
- // if (0 < f* < 10^(-x)) then the result is a midpoint
- // if floor(C*) is even then C* = floor(C*) - logical right
- // shift; C* has p decimal digits, correct by Prop. 1)
- // else if floor(C*) is odd C* = floor(C*)-1 (logical right
- // shift; C* has p decimal digits, correct by Pr. 1)
- // else
- // C* = floor(C*) (logical right shift; C has p decimal digits,
- // correct by Property 1)
- // n = C* * 10^(e+x)
- // shift right C* by Ex-64 = shiftright128[ind]
- shift = shiftright128[ind - 1]; // 0 <= shift <= 39
- Cstar = Cstar >> shift;
- // determine inexactness of the rounding of C*
- // if (0 < f* - 1/2 < 10^(-x)) then
- // the result is exact
- // else // if (f* - 1/2 > T*) then
- // the result is inexact
- if (ind - 1 <= 2) {
- if (fstar.w[0] > 0x8000000000000000ull) {
- // f* > 1/2 and the result may be exact
- tmp64 = fstar.w[0] - 0x8000000000000000ull; // f* - 1/2
- if ((tmp64 > ten2mk128trunc[ind - 1].w[1])) {
- // ten2mk128trunc[ind -1].w[1] is identical to
- // ten2mk128[ind -1].w[1]
- // set the inexact flag
- *pfpsf |= INEXACT_EXCEPTION;
- } // else the result is exact
- } else { // the result is inexact; f2* <= 1/2
- // set the inexact flag
- *pfpsf |= INEXACT_EXCEPTION;
- }
- } else { // if 3 <= ind - 1 <= 14
- if (fstar.w[1] > onehalf128[ind - 1] ||
- (fstar.w[1] == onehalf128[ind - 1] && fstar.w[0])) {
- // f2* > 1/2 and the result may be exact
- // Calculate f2* - 1/2
- tmp64 = fstar.w[1] - onehalf128[ind - 1];
- if (tmp64 || fstar.w[0] > ten2mk128trunc[ind - 1].w[1]) {
- // ten2mk128trunc[ind -1].w[1] is identical to
- // ten2mk128[ind -1].w[1]
- // set the inexact flag
- *pfpsf |= INEXACT_EXCEPTION;
- } // else the result is exact
- } else { // the result is inexact; f2* <= 1/2
- // set the inexact flag
- *pfpsf |= INEXACT_EXCEPTION;
- }
- }
- // if the result was a midpoint it was rounded away from zero, so
- // it will need a correction
- // check for midpoints
- if ((fstar.w[1] == 0) && fstar.w[0]
- && (fstar.w[0] <= ten2mk128trunc[ind - 1].w[1])) {
- // ten2mk128trunc[ind -1].w[1] is identical to
- // ten2mk128[ind -1].w[1]
- // the result is a midpoint; round to nearest
- if (Cstar & 0x01) { // Cstar is odd; MP in [EVEN, ODD]
- // if floor(C*) is odd C = floor(C*) - 1; the result >= 1
- Cstar--; // Cstar is now even
- } // else MP in [ODD, EVEN]
- }
- if (x_sign)
- res = -Cstar;
- else
- res = Cstar;
- } else if (exp == 0) {
- // 1 <= q <= 10
- // res = +/-C (exact)
- if (x_sign)
- res = -C1;
- else
- res = C1;
- } else { // if (exp > 0) => 1 <= exp <= 9, 1 <= q < 9, 2 <= q + exp <= 10
- // res = +/-C * 10^exp (exact)
- if (x_sign)
- res = -C1 * ten2k64[exp];
- else
- res = C1 * ten2k64[exp];
- }
- }
- BID_RETURN (res);
- }
- /*****************************************************************************
- * BID64_to_int32_floor
- ****************************************************************************/
- #if DECIMAL_CALL_BY_REFERENCE
- void
- bid64_to_int32_floor (int *pres,
- UINT64 *
- px _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
- _EXC_INFO_PARAM) {
- UINT64 x = *px;
- #else
- int
- bid64_to_int32_floor (UINT64 x _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
- _EXC_INFO_PARAM) {
- #endif
- int res;
- UINT64 x_sign;
- UINT64 x_exp;
- int exp; // unbiased exponent
- // Note: C1 represents x_significand (UINT64)
- UINT64 tmp64;
- BID_UI64DOUBLE tmp1;
- unsigned int x_nr_bits;
- int q, ind, shift;
- UINT64 C1;
- UINT64 Cstar; // C* represents up to 16 decimal digits ~ 54 bits
- UINT128 fstar;
- UINT128 P128;
- // check for NaN or Infinity
- if ((x & MASK_NAN) == MASK_NAN || (x & MASK_INF) == MASK_INF) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // unpack x
- x_sign = x & MASK_SIGN; // 0 for positive, MASK_SIGN for negative
- // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] =>
- if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) {
- x_exp = (x & MASK_BINARY_EXPONENT2) >> 51; // biased
- C1 = (x & MASK_BINARY_SIG2) | MASK_BINARY_OR2;
- if (C1 > 9999999999999999ull) { // non-canonical
- x_exp = 0;
- C1 = 0;
- }
- } else {
- x_exp = (x & MASK_BINARY_EXPONENT1) >> 53; // biased
- C1 = x & MASK_BINARY_SIG1;
- }
- // check for zeros (possibly from non-canonical values)
- if (C1 == 0x0ull) {
- // x is 0
- res = 0x00000000;
- BID_RETURN (res);
- }
- // x is not special and is not zero
- // q = nr. of decimal digits in x (1 <= q <= 54)
- // determine first the nr. of bits in x
- if (C1 >= 0x0020000000000000ull) { // x >= 2^53
- // split the 64-bit value in two 32-bit halves to avoid rounding errors
- if (C1 >= 0x0000000100000000ull) { // x >= 2^32
- tmp1.d = (double) (C1 >> 32); // exact conversion
- x_nr_bits =
- 33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- } else { // x < 2^32
- tmp1.d = (double) C1; // exact conversion
- x_nr_bits =
- 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- }
- } else { // if x < 2^53
- tmp1.d = (double) C1; // exact conversion
- x_nr_bits =
- 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- }
- q = nr_digits[x_nr_bits - 1].digits;
- if (q == 0) {
- q = nr_digits[x_nr_bits - 1].digits1;
- if (C1 >= nr_digits[x_nr_bits - 1].threshold_lo)
- q++;
- }
- exp = x_exp - 398; // unbiased exponent
- if ((q + exp) > 10) { // x >= 10^10 ~= 2^33.2... (cannot fit in 32 bits)
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- } else if ((q + exp) == 10) { // x = c(0)c(1)...c(9).c(10)...c(q-1)
- // in this case 2^29.89... ~= 10^9 <= x < 10^10 ~= 2^33.2...
- // so x rounded to an integer may or may not fit in a signed 32-bit int
- // the cases that do not fit are identified here; the ones that fit
- // fall through and will be handled with other cases further,
- // under '1 <= q + exp <= 10'
- if (x_sign) { // if n < 0 and q + exp = 10
- // if n < -2^31 then n is too large
- // too large if c(0)c(1)...c(9).c(10)...c(q-1) > 2^31
- // <=> 0.c(0)c(1)...c(q-1) * 10^11 > 0x500000000, 1<=q<=16
- // <=> C * 10^(11-q) >= 0x500000000, 1<=q<=16
- if (q <= 11) {
- // Note: C * 10^(11-q) has 10 or 11 digits; 0x500000000 has 11 digits
- tmp64 = C1 * ten2k64[11 - q]; // C scaled up to 11-digit int
- // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
- if (tmp64 > 0x500000000ull) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- } else { // if (q > 11), i.e. 12 <= q <= 16 and so -15 <= exp <= -2
- // C * 10^(11-q) > 0x500000000 <=>
- // C > 0x500000000 * 10^(q-11) where 1 <= q - 11 <= 5
- // (scale 2^31+1 up)
- // Note: 0x500000000*10^(q-11) has q-1 or q digits, where q <= 16
- tmp64 = 0x500000000ull * ten2k64[q - 11];
- if (C1 > tmp64) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- }
- } else { // if n > 0 and q + exp = 10
- // if n >= 2^31 then n is too large
- // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31
- // <=> 0.c(0)c(1)...c(q-1) * 10^11 >= 0x500000000, 1<=q<=16
- // <=> C * 10^(11-q) >= 0x500000000, 1<=q<=16
- if (q <= 11) {
- // Note: C * 10^(11-q) has 10 or 11 digits; 0x500000000 has 11 digits
- tmp64 = C1 * ten2k64[11 - q]; // C scaled up to 11-digit int
- // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
- if (tmp64 >= 0x500000000ull) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- } else { // if (q > 11), i.e. 12 <= q <= 16 and so -15 <= exp <= -2
- // C * 10^(11-q) >= 0x500000000 <=>
- // C >= 0x500000000 * 10^(q-11) where 1 <= q - 11 <= 5
- // (scale 2^31-1 up)
- // Note: 0x500000000*10^(q-11) has q-1 or q digits, where q <= 16
- tmp64 = 0x500000000ull * ten2k64[q - 11];
- if (C1 >= tmp64) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- }
- }
- }
- // n is not too large to be converted to int32: -2^31 <= n < 2^31
- // Note: some of the cases tested for above fall through to this point
- if ((q + exp) <= 0) { // n = +/-0.[0...0]c(0)c(1)...c(q-1)
- // return -1 or 0
- if (x_sign)
- res = 0xffffffff;
- else
- res = 0x00000000;
- BID_RETURN (res);
- } else { // if (1 <= q + exp <= 10, 1 <= q <= 16, -15 <= exp <= 9)
- // -2^31-1 < x <= -1 or 1 <= x < 2^31 so x can be rounded
- // to nearest to a 32-bit signed integer
- if (exp < 0) { // 2 <= q <= 16, -15 <= exp <= -1, 1 <= q + exp <= 10
- ind = -exp; // 1 <= ind <= 15; ind is a synonym for 'x'
- // chop off ind digits from the lower part of C1
- // C1 fits in 64 bits
- // calculate C* and f*
- // C* is actually floor(C*) in this case
- // C* and f* need shifting and masking, as shown by
- // shiftright128[] and maskhigh128[]
- // 1 <= x <= 15
- // kx = 10^(-x) = ten2mk64[ind - 1]
- // C* = C1 * 10^(-x)
- // the approximation of 10^(-x) was rounded up to 54 bits
- __mul_64x64_to_128MACH (P128, C1, ten2mk64[ind - 1]);
- Cstar = P128.w[1];
- fstar.w[1] = P128.w[1] & maskhigh128[ind - 1];
- fstar.w[0] = P128.w[0];
- // the top Ex bits of 10^(-x) are T* = ten2mk128trunc[ind].w[0], e.g.
- // if x=1, T*=ten2mk128trunc[0].w[0]=0x1999999999999999
- // C* = floor(C*) (logical right shift; C has p decimal digits,
- // correct by Property 1)
- // n = C* * 10^(e+x)
- // shift right C* by Ex-64 = shiftright128[ind]
- shift = shiftright128[ind - 1]; // 0 <= shift <= 39
- Cstar = Cstar >> shift;
- // determine inexactness of the rounding of C*
- // if (0 < f* < 10^(-x)) then
- // the result is exact
- // else // if (f* > T*) then
- // the result is inexact
- if (ind - 1 <= 2) {
- if (fstar.w[0] > ten2mk128trunc[ind - 1].w[1]) {
- // ten2mk128trunc[ind -1].w[1] is identical to
- // ten2mk128[ind -1].w[1]
- if (x_sign) { // negative and inexact
- Cstar++;
- }
- } // else the result is exact
- } else { // if 3 <= ind - 1 <= 14
- if (fstar.w[1] || fstar.w[0] > ten2mk128trunc[ind - 1].w[1]) {
- // ten2mk128trunc[ind -1].w[1] is identical to
- // ten2mk128[ind -1].w[1]
- if (x_sign) { // negative and inexact
- Cstar++;
- }
- } // else the result is exact
- }
- if (x_sign)
- res = -Cstar;
- else
- res = Cstar;
- } else if (exp == 0) {
- // 1 <= q <= 10
- // res = +/-C (exact)
- if (x_sign)
- res = -C1;
- else
- res = C1;
- } else { // if (exp > 0) => 1 <= exp <= 9, 1 <= q < 9, 2 <= q + exp <= 10
- // res = +/-C * 10^exp (exact)
- if (x_sign)
- res = -C1 * ten2k64[exp];
- else
- res = C1 * ten2k64[exp];
- }
- }
- BID_RETURN (res);
- }
- /*****************************************************************************
- * BID64_to_int32_xfloor
- ****************************************************************************/
- #if DECIMAL_CALL_BY_REFERENCE
- void
- bid64_to_int32_xfloor (int *pres,
- UINT64 *
- px _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
- _EXC_INFO_PARAM) {
- UINT64 x = *px;
- #else
- int
- bid64_to_int32_xfloor (UINT64 x _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
- _EXC_INFO_PARAM) {
- #endif
- int res;
- UINT64 x_sign;
- UINT64 x_exp;
- int exp; // unbiased exponent
- // Note: C1 represents x_significand (UINT64)
- UINT64 tmp64;
- BID_UI64DOUBLE tmp1;
- unsigned int x_nr_bits;
- int q, ind, shift;
- UINT64 C1;
- UINT64 Cstar; // C* represents up to 16 decimal digits ~ 54 bits
- UINT128 fstar;
- UINT128 P128;
- // check for NaN or Infinity
- if ((x & MASK_NAN) == MASK_NAN || (x & MASK_INF) == MASK_INF) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // unpack x
- x_sign = x & MASK_SIGN; // 0 for positive, MASK_SIGN for negative
- // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] =>
- if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) {
- x_exp = (x & MASK_BINARY_EXPONENT2) >> 51; // biased
- C1 = (x & MASK_BINARY_SIG2) | MASK_BINARY_OR2;
- if (C1 > 9999999999999999ull) { // non-canonical
- x_exp = 0;
- C1 = 0;
- }
- } else {
- x_exp = (x & MASK_BINARY_EXPONENT1) >> 53; // biased
- C1 = x & MASK_BINARY_SIG1;
- }
- // check for zeros (possibly from non-canonical values)
- if (C1 == 0x0ull) {
- // x is 0
- res = 0x00000000;
- BID_RETURN (res);
- }
- // x is not special and is not zero
- // q = nr. of decimal digits in x (1 <= q <= 54)
- // determine first the nr. of bits in x
- if (C1 >= 0x0020000000000000ull) { // x >= 2^53
- // split the 64-bit value in two 32-bit halves to avoid rounding errors
- if (C1 >= 0x0000000100000000ull) { // x >= 2^32
- tmp1.d = (double) (C1 >> 32); // exact conversion
- x_nr_bits =
- 33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- } else { // x < 2^32
- tmp1.d = (double) C1; // exact conversion
- x_nr_bits =
- 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- }
- } else { // if x < 2^53
- tmp1.d = (double) C1; // exact conversion
- x_nr_bits =
- 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- }
- q = nr_digits[x_nr_bits - 1].digits;
- if (q == 0) {
- q = nr_digits[x_nr_bits - 1].digits1;
- if (C1 >= nr_digits[x_nr_bits - 1].threshold_lo)
- q++;
- }
- exp = x_exp - 398; // unbiased exponent
- if ((q + exp) > 10) { // x >= 10^10 ~= 2^33.2... (cannot fit in 32 bits)
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- } else if ((q + exp) == 10) { // x = c(0)c(1)...c(9).c(10)...c(q-1)
- // in this case 2^29.89... ~= 10^9 <= x < 10^10 ~= 2^33.2...
- // so x rounded to an integer may or may not fit in a signed 32-bit int
- // the cases that do not fit are identified here; the ones that fit
- // fall through and will be handled with other cases further,
- // under '1 <= q + exp <= 10'
- if (x_sign) { // if n < 0 and q + exp = 10
- // if n < -2^31 then n is too large
- // too large if c(0)c(1)...c(9).c(10)...c(q-1) > 2^31
- // <=> 0.c(0)c(1)...c(q-1) * 10^11 > 0x500000000, 1<=q<=16
- // <=> C * 10^(11-q) >= 0x500000000, 1<=q<=16
- if (q <= 11) {
- // Note: C * 10^(11-q) has 10 or 11 digits; 0x500000000 has 11 digits
- tmp64 = C1 * ten2k64[11 - q]; // C scaled up to 11-digit int
- // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
- if (tmp64 > 0x500000000ull) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- } else { // if (q > 11), i.e. 12 <= q <= 16 and so -15 <= exp <= -2
- // C * 10^(11-q) > 0x500000000 <=>
- // C > 0x500000000 * 10^(q-11) where 1 <= q - 11 <= 5
- // (scale 2^31+1 up)
- // Note: 0x500000000*10^(q-11) has q-1 or q digits, where q <= 16
- tmp64 = 0x500000000ull * ten2k64[q - 11];
- if (C1 > tmp64) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- }
- } else { // if n > 0 and q + exp = 10
- // if n >= 2^31 then n is too large
- // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31
- // <=> 0.c(0)c(1)...c(q-1) * 10^11 >= 0x500000000, 1<=q<=16
- // <=> C * 10^(11-q) >= 0x500000000, 1<=q<=16
- if (q <= 11) {
- // Note: C * 10^(11-q) has 10 or 11 digits; 0x500000000 has 11 digits
- tmp64 = C1 * ten2k64[11 - q]; // C scaled up to 11-digit int
- // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
- if (tmp64 >= 0x500000000ull) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- } else { // if (q > 11), i.e. 12 <= q <= 16 and so -15 <= exp <= -2
- // C * 10^(11-q) >= 0x500000000 <=>
- // C >= 0x500000000 * 10^(q-11) where 1 <= q - 11 <= 5
- // (scale 2^31-1 up)
- // Note: 0x500000000*10^(q-11) has q-1 or q digits, where q <= 16
- tmp64 = 0x500000000ull * ten2k64[q - 11];
- if (C1 >= tmp64) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- }
- }
- }
- // n is not too large to be converted to int32: -2^31 <= n < 2^31
- // Note: some of the cases tested for above fall through to this point
- if ((q + exp) <= 0) { // n = +/-0.[0...0]c(0)c(1)...c(q-1)
- // set inexact flag
- *pfpsf |= INEXACT_EXCEPTION;
- // return -1 or 0
- if (x_sign)
- res = 0xffffffff;
- else
- res = 0x00000000;
- BID_RETURN (res);
- } else { // if (1 <= q + exp <= 10, 1 <= q <= 16, -15 <= exp <= 9)
- // -2^31-1 < x <= -1 or 1 <= x < 2^31 so x can be rounded
- // to nearest to a 32-bit signed integer
- if (exp < 0) { // 2 <= q <= 16, -15 <= exp <= -1, 1 <= q + exp <= 10
- ind = -exp; // 1 <= ind <= 15; ind is a synonym for 'x'
- // chop off ind digits from the lower part of C1
- // C1 fits in 64 bits
- // calculate C* and f*
- // C* is actually floor(C*) in this case
- // C* and f* need shifting and masking, as shown by
- // shiftright128[] and maskhigh128[]
- // 1 <= x <= 15
- // kx = 10^(-x) = ten2mk64[ind - 1]
- // C* = C1 * 10^(-x)
- // the approximation of 10^(-x) was rounded up to 54 bits
- __mul_64x64_to_128MACH (P128, C1, ten2mk64[ind - 1]);
- Cstar = P128.w[1];
- fstar.w[1] = P128.w[1] & maskhigh128[ind - 1];
- fstar.w[0] = P128.w[0];
- // the top Ex bits of 10^(-x) are T* = ten2mk128trunc[ind].w[0], e.g.
- // if x=1, T*=ten2mk128trunc[0].w[0]=0x1999999999999999
- // C* = floor(C*) (logical right shift; C has p decimal digits,
- // correct by Property 1)
- // n = C* * 10^(e+x)
- // shift right C* by Ex-64 = shiftright128[ind]
- shift = shiftright128[ind - 1]; // 0 <= shift <= 39
- Cstar = Cstar >> shift;
- // determine inexactness of the rounding of C*
- // if (0 < f* < 10^(-x)) then
- // the result is exact
- // else // if (f* > T*) then
- // the result is inexact
- if (ind - 1 <= 2) {
- if (fstar.w[0] > ten2mk128trunc[ind - 1].w[1]) {
- // ten2mk128trunc[ind -1].w[1] is identical to
- // ten2mk128[ind -1].w[1]
- if (x_sign) { // negative and inexact
- Cstar++;
- }
- // set the inexact flag
- *pfpsf |= INEXACT_EXCEPTION;
- } // else the result is exact
- } else { // if 3 <= ind - 1 <= 14
- if (fstar.w[1] || fstar.w[0] > ten2mk128trunc[ind - 1].w[1]) {
- // ten2mk128trunc[ind -1].w[1] is identical to
- // ten2mk128[ind -1].w[1]
- if (x_sign) { // negative and inexact
- Cstar++;
- }
- // set the inexact flag
- *pfpsf |= INEXACT_EXCEPTION;
- } // else the result is exact
- }
- if (x_sign)
- res = -Cstar;
- else
- res = Cstar;
- } else if (exp == 0) {
- // 1 <= q <= 10
- // res = +/-C (exact)
- if (x_sign)
- res = -C1;
- else
- res = C1;
- } else { // if (exp > 0) => 1 <= exp <= 9, 1 <= q < 9, 2 <= q + exp <= 10
- // res = +/-C * 10^exp (exact)
- if (x_sign)
- res = -C1 * ten2k64[exp];
- else
- res = C1 * ten2k64[exp];
- }
- }
- BID_RETURN (res);
- }
- /*****************************************************************************
- * BID64_to_int32_ceil
- ****************************************************************************/
- #if DECIMAL_CALL_BY_REFERENCE
- void
- bid64_to_int32_ceil (int *pres,
- UINT64 *
- px _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
- _EXC_INFO_PARAM) {
- UINT64 x = *px;
- #else
- int
- bid64_to_int32_ceil (UINT64 x _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
- _EXC_INFO_PARAM) {
- #endif
- int res;
- UINT64 x_sign;
- UINT64 x_exp;
- int exp; // unbiased exponent
- // Note: C1 represents x_significand (UINT64)
- UINT64 tmp64;
- BID_UI64DOUBLE tmp1;
- unsigned int x_nr_bits;
- int q, ind, shift;
- UINT64 C1;
- UINT64 Cstar; // C* represents up to 16 decimal digits ~ 54 bits
- UINT128 fstar;
- UINT128 P128;
- // check for NaN or Infinity
- if ((x & MASK_NAN) == MASK_NAN || (x & MASK_INF) == MASK_INF) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // unpack x
- x_sign = x & MASK_SIGN; // 0 for positive, MASK_SIGN for negative
- // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] =>
- if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) {
- x_exp = (x & MASK_BINARY_EXPONENT2) >> 51; // biased
- C1 = (x & MASK_BINARY_SIG2) | MASK_BINARY_OR2;
- if (C1 > 9999999999999999ull) { // non-canonical
- x_exp = 0;
- C1 = 0;
- }
- } else {
- x_exp = (x & MASK_BINARY_EXPONENT1) >> 53; // biased
- C1 = x & MASK_BINARY_SIG1;
- }
- // check for zeros (possibly from non-canonical values)
- if (C1 == 0x0ull) {
- // x is 0
- res = 0x00000000;
- BID_RETURN (res);
- }
- // x is not special and is not zero
- // q = nr. of decimal digits in x (1 <= q <= 54)
- // determine first the nr. of bits in x
- if (C1 >= 0x0020000000000000ull) { // x >= 2^53
- // split the 64-bit value in two 32-bit halves to avoid rounding errors
- if (C1 >= 0x0000000100000000ull) { // x >= 2^32
- tmp1.d = (double) (C1 >> 32); // exact conversion
- x_nr_bits =
- 33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- } else { // x < 2^32
- tmp1.d = (double) C1; // exact conversion
- x_nr_bits =
- 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- }
- } else { // if x < 2^53
- tmp1.d = (double) C1; // exact conversion
- x_nr_bits =
- 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- }
- q = nr_digits[x_nr_bits - 1].digits;
- if (q == 0) {
- q = nr_digits[x_nr_bits - 1].digits1;
- if (C1 >= nr_digits[x_nr_bits - 1].threshold_lo)
- q++;
- }
- exp = x_exp - 398; // unbiased exponent
- if ((q + exp) > 10) { // x >= 10^10 ~= 2^33.2... (cannot fit in 32 bits)
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- } else if ((q + exp) == 10) { // x = c(0)c(1)...c(9).c(10)...c(q-1)
- // in this case 2^29.89... ~= 10^9 <= x < 10^10 ~= 2^33.2...
- // so x rounded to an integer may or may not fit in a signed 32-bit int
- // the cases that do not fit are identified here; the ones that fit
- // fall through and will be handled with other cases further,
- // under '1 <= q + exp <= 10'
- if (x_sign) { // if n < 0 and q + exp = 10
- // if n <= -2^31 - 1 then n is too large
- // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31+1
- // <=> 0.c(0)c(1)...c(q-1) * 10^11 > 0x50000000a, 1<=q<=16
- // <=> C * 10^(11-q) >= 0x50000000a, 1<=q<=16
- if (q <= 11) {
- // Note: C * 10^(11-q) has 10 or 11 digits; 0x50000000a has 11 digits
- tmp64 = C1 * ten2k64[11 - q]; // C scaled up to 11-digit int
- // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
- if (tmp64 >= 0x50000000aull) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- } else { // if (q > 11), i.e. 12 <= q <= 16 and so -15 <= exp <= -2
- // C * 10^(11-q) >= 0x50000000a <=>
- // C >= 0x50000000a * 10^(q-11) where 1 <= q - 11 <= 5
- // (scale 2^31+1 up)
- // Note: 0x50000000a*10^(q-11) has q-1 or q digits, where q <= 16
- tmp64 = 0x50000000aull * ten2k64[q - 11];
- if (C1 >= tmp64) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- }
- } else { // if n > 0 and q + exp = 10
- // if n > 2^31 - 1 then n is too large
- // too large if c(0)c(1)...c(9).c(10)...c(q-1) > 2^31 - 1
- // <=> 0.c(0)c(1)...c(q-1) * 10^11 > 0x4fffffff6, 1<=q<=16
- // <=> C * 10^(11-q) > 0x4fffffff6, 1<=q<=16
- if (q <= 11) {
- // Note: C * 10^(11-q) has 10 or 11 digits; 0x4fffffff6 has 11 digits
- tmp64 = C1 * ten2k64[11 - q]; // C scaled up to 11-digit int
- // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
- if (tmp64 > 0x4fffffff6ull) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- } else { // if (q > 11), i.e. 12 <= q <= 16 and so -15 <= exp <= -2
- // C * 10^(11-q) > 0x4fffffff6 <=>
- // C > 0x4fffffff6 * 10^(q-11) where 1 <= q - 11 <= 5
- // (scale 2^31-1 up)
- // Note: 0x4fffffff6*10^(q-11) has q-1 or q digits, where q <= 16
- tmp64 = 0x4fffffff6ull * ten2k64[q - 11];
- if (C1 > tmp64) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- }
- }
- }
- // n is not too large to be converted to int32: -2^31 - 1 < n <= 2^31 - 1
- // Note: some of the cases tested for above fall through to this point
- if ((q + exp) <= 0) { // n = +/-0.[0...0]c(0)c(1)...c(q-1)
- // return 0 or 1
- if (x_sign)
- res = 0x00000000;
- else
- res = 0x00000001;
- BID_RETURN (res);
- } else { // if (1 <= q + exp <= 10, 1 <= q <= 16, -15 <= exp <= 9)
- // -2^31-1 < x <= -1 or 1 <= x <= 2^31-1 so x can be rounded
- // to nearest to a 32-bit signed integer
- if (exp < 0) { // 2 <= q <= 16, -15 <= exp <= -1, 1 <= q + exp <= 10
- ind = -exp; // 1 <= ind <= 15; ind is a synonym for 'x'
- // chop off ind digits from the lower part of C1
- // C1 fits in 64 bits
- // calculate C* and f*
- // C* is actually floor(C*) in this case
- // C* and f* need shifting and masking, as shown by
- // shiftright128[] and maskhigh128[]
- // 1 <= x <= 15
- // kx = 10^(-x) = ten2mk64[ind - 1]
- // C* = C1 * 10^(-x)
- // the approximation of 10^(-x) was rounded up to 54 bits
- __mul_64x64_to_128MACH (P128, C1, ten2mk64[ind - 1]);
- Cstar = P128.w[1];
- fstar.w[1] = P128.w[1] & maskhigh128[ind - 1];
- fstar.w[0] = P128.w[0];
- // the top Ex bits of 10^(-x) are T* = ten2mk128trunc[ind].w[0], e.g.
- // if x=1, T*=ten2mk128trunc[0].w[0]=0x1999999999999999
- // C* = floor(C*) (logical right shift; C has p decimal digits,
- // correct by Property 1)
- // n = C* * 10^(e+x)
- // shift right C* by Ex-64 = shiftright128[ind]
- shift = shiftright128[ind - 1]; // 0 <= shift <= 39
- Cstar = Cstar >> shift;
- // determine inexactness of the rounding of C*
- // if (0 < f* < 10^(-x)) then
- // the result is exact
- // else // if (f* > T*) then
- // the result is inexact
- if (ind - 1 <= 2) {
- if (fstar.w[0] > ten2mk128trunc[ind - 1].w[1]) {
- // ten2mk128trunc[ind -1].w[1] is identical to
- // ten2mk128[ind -1].w[1]
- if (!x_sign) { // positive and inexact
- Cstar++;
- }
- } // else the result is exact
- } else { // if 3 <= ind - 1 <= 14
- if (fstar.w[1] || fstar.w[0] > ten2mk128trunc[ind - 1].w[1]) {
- // ten2mk128trunc[ind -1].w[1] is identical to
- // ten2mk128[ind -1].w[1]
- if (!x_sign) { // positive and inexact
- Cstar++;
- }
- } // else the result is exact
- }
- if (x_sign)
- res = -Cstar;
- else
- res = Cstar;
- } else if (exp == 0) {
- // 1 <= q <= 10
- // res = +/-C (exact)
- if (x_sign)
- res = -C1;
- else
- res = C1;
- } else { // if (exp > 0) => 1 <= exp <= 9, 1 <= q < 9, 2 <= q + exp <= 10
- // res = +/-C * 10^exp (exact)
- if (x_sign)
- res = -C1 * ten2k64[exp];
- else
- res = C1 * ten2k64[exp];
- }
- }
- BID_RETURN (res);
- }
- /*****************************************************************************
- * BID64_to_int32_xceil
- ****************************************************************************/
- #if DECIMAL_CALL_BY_REFERENCE
- void
- bid64_to_int32_xceil (int *pres,
- UINT64 *
- px _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
- _EXC_INFO_PARAM) {
- UINT64 x = *px;
- #else
- int
- bid64_to_int32_xceil (UINT64 x _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
- _EXC_INFO_PARAM) {
- #endif
- int res;
- UINT64 x_sign;
- UINT64 x_exp;
- int exp; // unbiased exponent
- // Note: C1 represents x_significand (UINT64)
- UINT64 tmp64;
- BID_UI64DOUBLE tmp1;
- unsigned int x_nr_bits;
- int q, ind, shift;
- UINT64 C1;
- UINT64 Cstar; // C* represents up to 16 decimal digits ~ 54 bits
- UINT128 fstar;
- UINT128 P128;
- // check for NaN or Infinity
- if ((x & MASK_NAN) == MASK_NAN || (x & MASK_INF) == MASK_INF) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // unpack x
- x_sign = x & MASK_SIGN; // 0 for positive, MASK_SIGN for negative
- // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] =>
- if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) {
- x_exp = (x & MASK_BINARY_EXPONENT2) >> 51; // biased
- C1 = (x & MASK_BINARY_SIG2) | MASK_BINARY_OR2;
- if (C1 > 9999999999999999ull) { // non-canonical
- x_exp = 0;
- C1 = 0;
- }
- } else {
- x_exp = (x & MASK_BINARY_EXPONENT1) >> 53; // biased
- C1 = x & MASK_BINARY_SIG1;
- }
- // check for zeros (possibly from non-canonical values)
- if (C1 == 0x0ull) {
- // x is 0
- res = 0x00000000;
- BID_RETURN (res);
- }
- // x is not special and is not zero
- // q = nr. of decimal digits in x (1 <= q <= 54)
- // determine first the nr. of bits in x
- if (C1 >= 0x0020000000000000ull) { // x >= 2^53
- // split the 64-bit value in two 32-bit halves to avoid rounding errors
- if (C1 >= 0x0000000100000000ull) { // x >= 2^32
- tmp1.d = (double) (C1 >> 32); // exact conversion
- x_nr_bits =
- 33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- } else { // x < 2^32
- tmp1.d = (double) C1; // exact conversion
- x_nr_bits =
- 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- }
- } else { // if x < 2^53
- tmp1.d = (double) C1; // exact conversion
- x_nr_bits =
- 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- }
- q = nr_digits[x_nr_bits - 1].digits;
- if (q == 0) {
- q = nr_digits[x_nr_bits - 1].digits1;
- if (C1 >= nr_digits[x_nr_bits - 1].threshold_lo)
- q++;
- }
- exp = x_exp - 398; // unbiased exponent
- if ((q + exp) > 10) { // x >= 10^10 ~= 2^33.2... (cannot fit in 32 bits)
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- } else if ((q + exp) == 10) { // x = c(0)c(1)...c(9).c(10)...c(q-1)
- // in this case 2^29.89... ~= 10^9 <= x < 10^10 ~= 2^33.2...
- // so x rounded to an integer may or may not fit in a signed 32-bit int
- // the cases that do not fit are identified here; the ones that fit
- // fall through and will be handled with other cases further,
- // under '1 <= q + exp <= 10'
- if (x_sign) { // if n < 0 and q + exp = 10
- // if n <= -2^31 - 1 then n is too large
- // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31+1
- // <=> 0.c(0)c(1)...c(q-1) * 10^11 > 0x50000000a, 1<=q<=16
- // <=> C * 10^(11-q) >= 0x50000000a, 1<=q<=16
- if (q <= 11) {
- // Note: C * 10^(11-q) has 10 or 11 digits; 0x50000000a has 11 digits
- tmp64 = C1 * ten2k64[11 - q]; // C scaled up to 11-digit int
- // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
- if (tmp64 >= 0x50000000aull) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- } else { // if (q > 11), i.e. 12 <= q <= 16 and so -15 <= exp <= -2
- // C * 10^(11-q) >= 0x50000000a <=>
- // C >= 0x50000000a * 10^(q-11) where 1 <= q - 11 <= 5
- // (scale 2^31+1 up)
- // Note: 0x50000000a*10^(q-11) has q-1 or q digits, where q <= 16
- tmp64 = 0x50000000aull * ten2k64[q - 11];
- if (C1 >= tmp64) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- }
- } else { // if n > 0 and q + exp = 10
- // if n > 2^31 - 1 then n is too large
- // too large if c(0)c(1)...c(9).c(10)...c(q-1) > 2^31 - 1
- // <=> 0.c(0)c(1)...c(q-1) * 10^11 > 0x4fffffff6, 1<=q<=16
- // <=> C * 10^(11-q) > 0x4fffffff6, 1<=q<=16
- if (q <= 11) {
- // Note: C * 10^(11-q) has 10 or 11 digits; 0x4fffffff6 has 11 digits
- tmp64 = C1 * ten2k64[11 - q]; // C scaled up to 11-digit int
- // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
- if (tmp64 > 0x4fffffff6ull) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- } else { // if (q > 11), i.e. 12 <= q <= 16 and so -15 <= exp <= -2
- // C * 10^(11-q) > 0x4fffffff6 <=>
- // C > 0x4fffffff6 * 10^(q-11) where 1 <= q - 11 <= 5
- // (scale 2^31-1 up)
- // Note: 0x4fffffff6*10^(q-11) has q-1 or q digits, where q <= 16
- tmp64 = 0x4fffffff6ull * ten2k64[q - 11];
- if (C1 > tmp64) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- }
- }
- }
- // n is not too large to be converted to int32: -2^31 - 1 < n <= 2^31 - 1
- // Note: some of the cases tested for above fall through to this point
- if ((q + exp) <= 0) { // n = +/-0.[0...0]c(0)c(1)...c(q-1)
- // set inexact flag
- *pfpsf |= INEXACT_EXCEPTION;
- // return 0 or 1
- if (x_sign)
- res = 0x00000000;
- else
- res = 0x00000001;
- BID_RETURN (res);
- } else { // if (1 <= q + exp <= 10, 1 <= q <= 16, -15 <= exp <= 9)
- // -2^31-1 < x <= -1 or 1 <= x <= 2^31-1 so x can be rounded
- // to nearest to a 32-bit signed integer
- if (exp < 0) { // 2 <= q <= 16, -15 <= exp <= -1, 1 <= q + exp <= 10
- ind = -exp; // 1 <= ind <= 15; ind is a synonym for 'x'
- // chop off ind digits from the lower part of C1
- // C1 fits in 64 bits
- // calculate C* and f*
- // C* is actually floor(C*) in this case
- // C* and f* need shifting and masking, as shown by
- // shiftright128[] and maskhigh128[]
- // 1 <= x <= 15
- // kx = 10^(-x) = ten2mk64[ind - 1]
- // C* = C1 * 10^(-x)
- // the approximation of 10^(-x) was rounded up to 54 bits
- __mul_64x64_to_128MACH (P128, C1, ten2mk64[ind - 1]);
- Cstar = P128.w[1];
- fstar.w[1] = P128.w[1] & maskhigh128[ind - 1];
- fstar.w[0] = P128.w[0];
- // the top Ex bits of 10^(-x) are T* = ten2mk128trunc[ind].w[0], e.g.
- // if x=1, T*=ten2mk128trunc[0].w[0]=0x1999999999999999
- // C* = floor(C*) (logical right shift; C has p decimal digits,
- // correct by Property 1)
- // n = C* * 10^(e+x)
- // shift right C* by Ex-64 = shiftright128[ind]
- shift = shiftright128[ind - 1]; // 0 <= shift <= 39
- Cstar = Cstar >> shift;
- // determine inexactness of the rounding of C*
- // if (0 < f* < 10^(-x)) then
- // the result is exact
- // else // if (f* > T*) then
- // the result is inexact
- if (ind - 1 <= 2) {
- if (fstar.w[0] > ten2mk128trunc[ind - 1].w[1]) {
- // ten2mk128trunc[ind -1].w[1] is identical to
- // ten2mk128[ind -1].w[1]
- if (!x_sign) { // positive and inexact
- Cstar++;
- }
- // set the inexact flag
- *pfpsf |= INEXACT_EXCEPTION;
- } // else the result is exact
- } else { // if 3 <= ind - 1 <= 14
- if (fstar.w[1] || fstar.w[0] > ten2mk128trunc[ind - 1].w[1]) {
- // ten2mk128trunc[ind -1].w[1] is identical to
- // ten2mk128[ind -1].w[1]
- if (!x_sign) { // positive and inexact
- Cstar++;
- }
- // set the inexact flag
- *pfpsf |= INEXACT_EXCEPTION;
- } // else the result is exact
- }
- if (x_sign)
- res = -Cstar;
- else
- res = Cstar;
- } else if (exp == 0) {
- // 1 <= q <= 10
- // res = +/-C (exact)
- if (x_sign)
- res = -C1;
- else
- res = C1;
- } else { // if (exp > 0) => 1 <= exp <= 9, 1 <= q < 9, 2 <= q + exp <= 10
- // res = +/-C * 10^exp (exact)
- if (x_sign)
- res = -C1 * ten2k64[exp];
- else
- res = C1 * ten2k64[exp];
- }
- }
- BID_RETURN (res);
- }
- /*****************************************************************************
- * BID64_to_int32_int
- ****************************************************************************/
- #if DECIMAL_CALL_BY_REFERENCE
- void
- bid64_to_int32_int (int *pres,
- UINT64 *
- px _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
- _EXC_INFO_PARAM) {
- UINT64 x = *px;
- #else
- int
- bid64_to_int32_int (UINT64 x _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
- _EXC_INFO_PARAM) {
- #endif
- int res;
- UINT64 x_sign;
- UINT64 x_exp;
- int exp; // unbiased exponent
- // Note: C1 represents x_significand (UINT64)
- UINT64 tmp64;
- BID_UI64DOUBLE tmp1;
- unsigned int x_nr_bits;
- int q, ind, shift;
- UINT64 C1;
- UINT64 Cstar; // C* represents up to 16 decimal digits ~ 54 bits
- UINT128 P128;
- // check for NaN or Infinity
- if ((x & MASK_NAN) == MASK_NAN || (x & MASK_INF) == MASK_INF) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // unpack x
- x_sign = x & MASK_SIGN; // 0 for positive, MASK_SIGN for negative
- // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] =>
- if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) {
- x_exp = (x & MASK_BINARY_EXPONENT2) >> 51; // biased
- C1 = (x & MASK_BINARY_SIG2) | MASK_BINARY_OR2;
- if (C1 > 9999999999999999ull) { // non-canonical
- x_exp = 0;
- C1 = 0;
- }
- } else {
- x_exp = (x & MASK_BINARY_EXPONENT1) >> 53; // biased
- C1 = x & MASK_BINARY_SIG1;
- }
- // check for zeros (possibly from non-canonical values)
- if (C1 == 0x0ull) {
- // x is 0
- res = 0x00000000;
- BID_RETURN (res);
- }
- // x is not special and is not zero
- // q = nr. of decimal digits in x (1 <= q <= 54)
- // determine first the nr. of bits in x
- if (C1 >= 0x0020000000000000ull) { // x >= 2^53
- // split the 64-bit value in two 32-bit halves to avoid rounding errors
- if (C1 >= 0x0000000100000000ull) { // x >= 2^32
- tmp1.d = (double) (C1 >> 32); // exact conversion
- x_nr_bits =
- 33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- } else { // x < 2^32
- tmp1.d = (double) C1; // exact conversion
- x_nr_bits =
- 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- }
- } else { // if x < 2^53
- tmp1.d = (double) C1; // exact conversion
- x_nr_bits =
- 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- }
- q = nr_digits[x_nr_bits - 1].digits;
- if (q == 0) {
- q = nr_digits[x_nr_bits - 1].digits1;
- if (C1 >= nr_digits[x_nr_bits - 1].threshold_lo)
- q++;
- }
- exp = x_exp - 398; // unbiased exponent
- if ((q + exp) > 10) { // x >= 10^10 ~= 2^33.2... (cannot fit in 32 bits)
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- } else if ((q + exp) == 10) { // x = c(0)c(1)...c(9).c(10)...c(q-1)
- // in this case 2^29.89... ~= 10^9 <= x < 10^10 ~= 2^33.2...
- // so x rounded to an integer may or may not fit in a signed 32-bit int
- // the cases that do not fit are identified here; the ones that fit
- // fall through and will be handled with other cases further,
- // under '1 <= q + exp <= 10'
- if (x_sign) { // if n < 0 and q + exp = 10
- // if n <= -2^31 - 1 then n is too large
- // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31+1
- // <=> 0.c(0)c(1)...c(q-1) * 10^11 > 0x50000000a, 1<=q<=16
- // <=> C * 10^(11-q) >= 0x50000000a, 1<=q<=16
- if (q <= 11) {
- // Note: C * 10^(11-q) has 10 or 11 digits; 0x50000000a has 11 digits
- tmp64 = C1 * ten2k64[11 - q]; // C scaled up to 11-digit int
- // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
- if (tmp64 >= 0x50000000aull) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- } else { // if (q > 11), i.e. 12 <= q <= 16 and so -15 <= exp <= -2
- // C * 10^(11-q) >= 0x50000000a <=>
- // C >= 0x50000000a * 10^(q-11) where 1 <= q - 11 <= 5
- // (scale 2^31+1 up)
- // Note: 0x50000000a*10^(q-11) has q-1 or q digits, where q <= 16
- tmp64 = 0x50000000aull * ten2k64[q - 11];
- if (C1 >= tmp64) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- }
- } else { // if n > 0 and q + exp = 10
- // if n >= 2^31 then n is too large
- // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31
- // <=> 0.c(0)c(1)...c(q-1) * 10^11 >= 0x500000000, 1<=q<=16
- // <=> C * 10^(11-q) >= 0x500000000, 1<=q<=16
- if (q <= 11) {
- // Note: C * 10^(11-q) has 10 or 11 digits; 0x500000000 has 11 digits
- tmp64 = C1 * ten2k64[11 - q]; // C scaled up to 11-digit int
- // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
- if (tmp64 >= 0x500000000ull) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- } else { // if (q > 11), i.e. 12 <= q <= 16 and so -15 <= exp <= -2
- // C * 10^(11-q) >= 0x500000000 <=>
- // C >= 0x500000000 * 10^(q-11) where 1 <= q - 11 <= 5
- // (scale 2^31-1 up)
- // Note: 0x500000000*10^(q-11) has q-1 or q digits, where q <= 16
- tmp64 = 0x500000000ull * ten2k64[q - 11];
- if (C1 >= tmp64) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- }
- }
- }
- // n is not too large to be converted to int32: -2^31 - 1 < n < 2^31
- // Note: some of the cases tested for above fall through to this point
- if ((q + exp) <= 0) { // n = +/-0.[0...0]c(0)c(1)...c(q-1)
- // return 0
- res = 0x00000000;
- BID_RETURN (res);
- } else { // if (1 <= q + exp <= 10, 1 <= q <= 16, -15 <= exp <= 9)
- // -2^31-1 < x <= -1 or 1 <= x < 2^31 so x can be rounded
- // to nearest to a 32-bit signed integer
- if (exp < 0) { // 2 <= q <= 16, -15 <= exp <= -1, 1 <= q + exp <= 10
- ind = -exp; // 1 <= ind <= 15; ind is a synonym for 'x'
- // chop off ind digits from the lower part of C1
- // C1 fits in 64 bits
- // calculate C* and f*
- // C* is actually floor(C*) in this case
- // C* and f* need shifting and masking, as shown by
- // shiftright128[] and maskhigh128[]
- // 1 <= x <= 15
- // kx = 10^(-x) = ten2mk64[ind - 1]
- // C* = C1 * 10^(-x)
- // the approximation of 10^(-x) was rounded up to 54 bits
- __mul_64x64_to_128MACH (P128, C1, ten2mk64[ind - 1]);
- Cstar = P128.w[1];
- // the top Ex bits of 10^(-x) are T* = ten2mk128trunc[ind].w[0], e.g.
- // if x=1, T*=ten2mk128trunc[0].w[0]=0x1999999999999999
- // C* = floor(C*) (logical right shift; C has p decimal digits,
- // correct by Property 1)
- // n = C* * 10^(e+x)
- // shift right C* by Ex-64 = shiftright128[ind]
- shift = shiftright128[ind - 1]; // 0 <= shift <= 39
- Cstar = Cstar >> shift;
- if (x_sign)
- res = -Cstar;
- else
- res = Cstar;
- } else if (exp == 0) {
- // 1 <= q <= 10
- // res = +/-C (exact)
- if (x_sign)
- res = -C1;
- else
- res = C1;
- } else { // if (exp > 0) => 1 <= exp <= 9, 1 <= q < 9, 2 <= q + exp <= 10
- // res = +/-C * 10^exp (exact)
- if (x_sign)
- res = -C1 * ten2k64[exp];
- else
- res = C1 * ten2k64[exp];
- }
- }
- BID_RETURN (res);
- }
- /*****************************************************************************
- * BID64_to_int32_xint
- ****************************************************************************/
- #if DECIMAL_CALL_BY_REFERENCE
- void
- bid64_to_int32_xint (int *pres,
- UINT64 *
- px _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
- _EXC_INFO_PARAM) {
- UINT64 x = *px;
- #else
- int
- bid64_to_int32_xint (UINT64 x _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
- _EXC_INFO_PARAM) {
- #endif
- int res;
- UINT64 x_sign;
- UINT64 x_exp;
- int exp; // unbiased exponent
- // Note: C1 represents x_significand (UINT64)
- UINT64 tmp64;
- BID_UI64DOUBLE tmp1;
- unsigned int x_nr_bits;
- int q, ind, shift;
- UINT64 C1;
- UINT64 Cstar; // C* represents up to 16 decimal digits ~ 54 bits
- UINT128 fstar;
- UINT128 P128;
- // check for NaN or Infinity
- if ((x & MASK_NAN) == MASK_NAN || (x & MASK_INF) == MASK_INF) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // unpack x
- x_sign = x & MASK_SIGN; // 0 for positive, MASK_SIGN for negative
- // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] =>
- if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) {
- x_exp = (x & MASK_BINARY_EXPONENT2) >> 51; // biased
- C1 = (x & MASK_BINARY_SIG2) | MASK_BINARY_OR2;
- if (C1 > 9999999999999999ull) { // non-canonical
- x_exp = 0;
- C1 = 0;
- }
- } else {
- x_exp = (x & MASK_BINARY_EXPONENT1) >> 53; // biased
- C1 = x & MASK_BINARY_SIG1;
- }
- // check for zeros (possibly from non-canonical values)
- if (C1 == 0x0ull) {
- // x is 0
- res = 0x00000000;
- BID_RETURN (res);
- }
- // x is not special and is not zero
- // q = nr. of decimal digits in x (1 <= q <= 54)
- // determine first the nr. of bits in x
- if (C1 >= 0x0020000000000000ull) { // x >= 2^53
- // split the 64-bit value in two 32-bit halves to avoid rounding errors
- if (C1 >= 0x0000000100000000ull) { // x >= 2^32
- tmp1.d = (double) (C1 >> 32); // exact conversion
- x_nr_bits =
- 33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- } else { // x < 2^32
- tmp1.d = (double) C1; // exact conversion
- x_nr_bits =
- 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- }
- } else { // if x < 2^53
- tmp1.d = (double) C1; // exact conversion
- x_nr_bits =
- 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- }
- q = nr_digits[x_nr_bits - 1].digits;
- if (q == 0) {
- q = nr_digits[x_nr_bits - 1].digits1;
- if (C1 >= nr_digits[x_nr_bits - 1].threshold_lo)
- q++;
- }
- exp = x_exp - 398; // unbiased exponent
- if ((q + exp) > 10) { // x >= 10^10 ~= 2^33.2... (cannot fit in 32 bits)
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- } else if ((q + exp) == 10) { // x = c(0)c(1)...c(9).c(10)...c(q-1)
- // in this case 2^29.89... ~= 10^9 <= x < 10^10 ~= 2^33.2...
- // so x rounded to an integer may or may not fit in a signed 32-bit int
- // the cases that do not fit are identified here; the ones that fit
- // fall through and will be handled with other cases further,
- // under '1 <= q + exp <= 10'
- if (x_sign) { // if n < 0 and q + exp = 10
- // if n <= -2^31 - 1 then n is too large
- // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31+1
- // <=> 0.c(0)c(1)...c(q-1) * 10^11 > 0x50000000a, 1<=q<=16
- // <=> C * 10^(11-q) >= 0x50000000a, 1<=q<=16
- if (q <= 11) {
- // Note: C * 10^(11-q) has 10 or 11 digits; 0x50000000a has 11 digits
- tmp64 = C1 * ten2k64[11 - q]; // C scaled up to 11-digit int
- // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
- if (tmp64 >= 0x50000000aull) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- } else { // if (q > 11), i.e. 12 <= q <= 16 and so -15 <= exp <= -2
- // C * 10^(11-q) >= 0x50000000a <=>
- // C >= 0x50000000a * 10^(q-11) where 1 <= q - 11 <= 5
- // (scale 2^31+1 up)
- // Note: 0x50000000a*10^(q-11) has q-1 or q digits, where q <= 16
- tmp64 = 0x50000000aull * ten2k64[q - 11];
- if (C1 >= tmp64) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- }
- } else { // if n > 0 and q + exp = 10
- // if n >= 2^31 then n is too large
- // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31
- // <=> 0.c(0)c(1)...c(q-1) * 10^11 >= 0x500000000, 1<=q<=16
- // <=> C * 10^(11-q) >= 0x500000000, 1<=q<=16
- if (q <= 11) {
- // Note: C * 10^(11-q) has 10 or 11 digits; 0x500000000 has 11 digits
- tmp64 = C1 * ten2k64[11 - q]; // C scaled up to 11-digit int
- // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
- if (tmp64 >= 0x500000000ull) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- } else { // if (q > 11), i.e. 12 <= q <= 16 and so -15 <= exp <= -2
- // C * 10^(11-q) >= 0x500000000 <=>
- // C >= 0x500000000 * 10^(q-11) where 1 <= q - 11 <= 5
- // (scale 2^31-1 up)
- // Note: 0x500000000*10^(q-11) has q-1 or q digits, where q <= 16
- tmp64 = 0x500000000ull * ten2k64[q - 11];
- if (C1 >= tmp64) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- }
- }
- }
- // n is not too large to be converted to int32: -2^31 - 1 < n < 2^31
- // Note: some of the cases tested for above fall through to this point
- if ((q + exp) <= 0) { // n = +/-0.[0...0]c(0)c(1)...c(q-1)
- // set inexact flag
- *pfpsf |= INEXACT_EXCEPTION;
- // return 0
- res = 0x00000000;
- BID_RETURN (res);
- } else { // if (1 <= q + exp <= 10, 1 <= q <= 16, -15 <= exp <= 9)
- // -2^31-1 < x <= -1 or 1 <= x < 2^31 so x can be rounded
- // to nearest to a 32-bit signed integer
- if (exp < 0) { // 2 <= q <= 16, -15 <= exp <= -1, 1 <= q + exp <= 10
- ind = -exp; // 1 <= ind <= 15; ind is a synonym for 'x'
- // chop off ind digits from the lower part of C1
- // C1 fits in 64 bits
- // calculate C* and f*
- // C* is actually floor(C*) in this case
- // C* and f* need shifting and masking, as shown by
- // shiftright128[] and maskhigh128[]
- // 1 <= x <= 15
- // kx = 10^(-x) = ten2mk64[ind - 1]
- // C* = C1 * 10^(-x)
- // the approximation of 10^(-x) was rounded up to 54 bits
- __mul_64x64_to_128MACH (P128, C1, ten2mk64[ind - 1]);
- Cstar = P128.w[1];
- fstar.w[1] = P128.w[1] & maskhigh128[ind - 1];
- fstar.w[0] = P128.w[0];
- // the top Ex bits of 10^(-x) are T* = ten2mk128trunc[ind].w[0], e.g.
- // if x=1, T*=ten2mk128trunc[0].w[0]=0x1999999999999999
- // C* = floor(C*) (logical right shift; C has p decimal digits,
- // correct by Property 1)
- // n = C* * 10^(e+x)
- // shift right C* by Ex-64 = shiftright128[ind]
- shift = shiftright128[ind - 1]; // 0 <= shift <= 39
- Cstar = Cstar >> shift;
- // determine inexactness of the rounding of C*
- // if (0 < f* < 10^(-x)) then
- // the result is exact
- // else // if (f* > T*) then
- // the result is inexact
- if (ind - 1 <= 2) {
- if (fstar.w[0] > ten2mk128trunc[ind - 1].w[1]) {
- // ten2mk128trunc[ind -1].w[1] is identical to
- // ten2mk128[ind -1].w[1]
- // set the inexact flag
- *pfpsf |= INEXACT_EXCEPTION;
- } // else the result is exact
- } else { // if 3 <= ind - 1 <= 14
- if (fstar.w[1] || fstar.w[0] > ten2mk128trunc[ind - 1].w[1]) {
- // ten2mk128trunc[ind -1].w[1] is identical to
- // ten2mk128[ind -1].w[1]
- // set the inexact flag
- *pfpsf |= INEXACT_EXCEPTION;
- } // else the result is exact
- }
- if (x_sign)
- res = -Cstar;
- else
- res = Cstar;
- } else if (exp == 0) {
- // 1 <= q <= 10
- // res = +/-C (exact)
- if (x_sign)
- res = -C1;
- else
- res = C1;
- } else { // if (exp > 0) => 1 <= exp <= 9, 1 <= q < 9, 2 <= q + exp <= 10
- // res = +/-C * 10^exp (exact)
- if (x_sign)
- res = -C1 * ten2k64[exp];
- else
- res = C1 * ten2k64[exp];
- }
- }
- BID_RETURN (res);
- }
- /*****************************************************************************
- * BID64_to_int32_rninta
- ****************************************************************************/
- #if DECIMAL_CALL_BY_REFERENCE
- void
- bid64_to_int32_rninta (int *pres,
- UINT64 *
- px _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
- _EXC_INFO_PARAM) {
- UINT64 x = *px;
- #else
- int
- bid64_to_int32_rninta (UINT64 x _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
- _EXC_INFO_PARAM) {
- #endif
- int res;
- UINT64 x_sign;
- UINT64 x_exp;
- int exp; // unbiased exponent
- // Note: C1 represents x_significand (UINT64)
- UINT64 tmp64;
- BID_UI64DOUBLE tmp1;
- unsigned int x_nr_bits;
- int q, ind, shift;
- UINT64 C1;
- UINT64 Cstar; // C* represents up to 16 decimal digits ~ 54 bits
- UINT128 P128;
- // check for NaN or Infinity
- if ((x & MASK_NAN) == MASK_NAN || (x & MASK_INF) == MASK_INF) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // unpack x
- x_sign = x & MASK_SIGN; // 0 for positive, MASK_SIGN for negative
- // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] =>
- if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) {
- x_exp = (x & MASK_BINARY_EXPONENT2) >> 51; // biased
- C1 = (x & MASK_BINARY_SIG2) | MASK_BINARY_OR2;
- if (C1 > 9999999999999999ull) { // non-canonical
- x_exp = 0;
- C1 = 0;
- }
- } else {
- x_exp = (x & MASK_BINARY_EXPONENT1) >> 53; // biased
- C1 = x & MASK_BINARY_SIG1;
- }
- // check for zeros (possibly from non-canonical values)
- if (C1 == 0x0ull) {
- // x is 0
- res = 0x00000000;
- BID_RETURN (res);
- }
- // x is not special and is not zero
- // q = nr. of decimal digits in x (1 <= q <= 54)
- // determine first the nr. of bits in x
- if (C1 >= 0x0020000000000000ull) { // x >= 2^53
- // split the 64-bit value in two 32-bit halves to avoid rounding errors
- if (C1 >= 0x0000000100000000ull) { // x >= 2^32
- tmp1.d = (double) (C1 >> 32); // exact conversion
- x_nr_bits =
- 33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- } else { // x < 2^32
- tmp1.d = (double) C1; // exact conversion
- x_nr_bits =
- 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- }
- } else { // if x < 2^53
- tmp1.d = (double) C1; // exact conversion
- x_nr_bits =
- 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- }
- q = nr_digits[x_nr_bits - 1].digits;
- if (q == 0) {
- q = nr_digits[x_nr_bits - 1].digits1;
- if (C1 >= nr_digits[x_nr_bits - 1].threshold_lo)
- q++;
- }
- exp = x_exp - 398; // unbiased exponent
- if ((q + exp) > 10) { // x >= 10^10 ~= 2^33.2... (cannot fit in 32 bits)
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- } else if ((q + exp) == 10) { // x = c(0)c(1)...c(9).c(10)...c(q-1)
- // in this case 2^29.89... ~= 10^9 <= x < 10^10 ~= 2^33.2...
- // so x rounded to an integer may or may not fit in a signed 32-bit int
- // the cases that do not fit are identified here; the ones that fit
- // fall through and will be handled with other cases further,
- // under '1 <= q + exp <= 10'
- if (x_sign) { // if n < 0 and q + exp = 10
- // if n <= -2^31 - 1/2 then n is too large
- // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31+1/2
- // <=> 0.c(0)c(1)...c(q-1) * 10^11 >= 0x500000005, 1<=q<=16
- // <=> C * 10^(11-q) >= 0x500000005, 1<=q<=16
- if (q <= 11) {
- // Note: C * 10^(11-q) has 10 or 11 digits; 0x500000005 has 11 digits
- tmp64 = C1 * ten2k64[11 - q]; // C scaled up to 11-digit int
- // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
- if (tmp64 >= 0x500000005ull) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- } else { // if (q > 11), i.e. 12 <= q <= 16 and so -15 <= exp <= -2
- // C * 10^(11-q) >= 0x500000005 <=>
- // C >= 0x500000005 * 10^(q-11) where 1 <= q - 11 <= 5
- // (scale 2^31+1/2 up)
- // Note: 0x500000005*10^(q-11) has q-1 or q digits, where q <= 16
- tmp64 = 0x500000005ull * ten2k64[q - 11];
- if (C1 >= tmp64) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- }
- } else { // if n > 0 and q + exp = 10
- // if n >= 2^31 - 1/2 then n is too large
- // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31-1/2
- // <=> 0.c(0)c(1)...c(q-1) * 10^11 >= 0x4fffffffb, 1<=q<=16
- // <=> C * 10^(11-q) >= 0x4fffffffb, 1<=q<=16
- if (q <= 11) {
- // Note: C * 10^(11-q) has 10 or 11 digits; 0x500000005 has 11 digits
- tmp64 = C1 * ten2k64[11 - q]; // C scaled up to 11-digit int
- // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
- if (tmp64 >= 0x4fffffffbull) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- } else { // if (q > 11), i.e. 12 <= q <= 16 and so -15 <= exp <= -2
- // C * 10^(11-q) >= 0x4fffffffb <=>
- // C >= 0x4fffffffb * 10^(q-11) where 1 <= q - 11 <= 5
- // (scale 2^31-1/2 up)
- // Note: 0x4fffffffb*10^(q-11) has q-1 or q digits, where q <= 16
- tmp64 = 0x4fffffffbull * ten2k64[q - 11];
- if (C1 >= tmp64) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- }
- }
- }
- // n is not too large to be converted to int32: -2^31 - 1/2 < n < 2^31 - 1/2
- // Note: some of the cases tested for above fall through to this point
- if ((q + exp) < 0) { // n = +/-0.0...c(0)c(1)...c(q-1)
- // return 0
- res = 0x00000000;
- BID_RETURN (res);
- } else if ((q + exp) == 0) { // n = +/-0.c(0)c(1)...c(q-1)
- // if 0.c(0)c(1)...c(q-1) <= 0.5 <=> c(0)c(1)...c(q-1) <= 5 * 10^(q-1)
- // res = 0
- // else
- // res = +/-1
- ind = q - 1;
- if (C1 < midpoint64[ind]) {
- res = 0x00000000; // return 0
- } else if (x_sign) { // n < 0
- res = 0xffffffff; // return -1
- } else { // n > 0
- res = 0x00000001; // return +1
- }
- } else { // if (1 <= q + exp <= 10, 1 <= q <= 16, -15 <= exp <= 9)
- // -2^31-1/2 <= x <= -1 or 1 <= x < 2^31-1/2 so x can be rounded
- // to nearest away to a 32-bit signed integer
- if (exp < 0) { // 2 <= q <= 16, -15 <= exp <= -1, 1 <= q + exp <= 10
- ind = -exp; // 1 <= ind <= 15; ind is a synonym for 'x'
- // chop off ind digits from the lower part of C1
- // C1 = C1 + 1/2 * 10^ind where the result C1 fits in 64 bits
- C1 = C1 + midpoint64[ind - 1];
- // calculate C* and f*
- // C* is actually floor(C*) in this case
- // C* and f* need shifting and masking, as shown by
- // shiftright128[] and maskhigh128[]
- // 1 <= x <= 15
- // kx = 10^(-x) = ten2mk64[ind - 1]
- // C* = (C1 + 1/2 * 10^x) * 10^(-x)
- // the approximation of 10^(-x) was rounded up to 54 bits
- __mul_64x64_to_128MACH (P128, C1, ten2mk64[ind - 1]);
- Cstar = P128.w[1];
- // the top Ex bits of 10^(-x) are T* = ten2mk128trunc[ind].w[0], e.g.
- // if x=1, T*=ten2mk128trunc[0].w[0]=0x1999999999999999
- // C* = floor(C*)-1 (logical right shift; C* has p decimal digits,
- // correct by Pr. 1)
- // n = C* * 10^(e+x)
- // shift right C* by Ex-64 = shiftright128[ind]
- shift = shiftright128[ind - 1]; // 0 <= shift <= 39
- Cstar = Cstar >> shift;
- // if the result was a midpoint it was rounded away from zero
- if (x_sign)
- res = -Cstar;
- else
- res = Cstar;
- } else if (exp == 0) {
- // 1 <= q <= 10
- // res = +/-C (exact)
- if (x_sign)
- res = -C1;
- else
- res = C1;
- } else { // if (exp > 0) => 1 <= exp <= 9, 1 <= q < 9, 2 <= q + exp <= 10
- // res = +/-C * 10^exp (exact)
- if (x_sign)
- res = -C1 * ten2k64[exp];
- else
- res = C1 * ten2k64[exp];
- }
- }
- BID_RETURN (res);
- }
- /*****************************************************************************
- * BID64_to_int32_xrninta
- ****************************************************************************/
- #if DECIMAL_CALL_BY_REFERENCE
- void
- bid64_to_int32_xrninta (int *pres,
- UINT64 *
- px _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
- _EXC_INFO_PARAM) {
- UINT64 x = *px;
- #else
- int
- bid64_to_int32_xrninta (UINT64 x _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
- _EXC_INFO_PARAM) {
- #endif
- int res;
- UINT64 x_sign;
- UINT64 x_exp;
- int exp; // unbiased exponent
- // Note: C1 represents x_significand (UINT64)
- UINT64 tmp64;
- BID_UI64DOUBLE tmp1;
- unsigned int x_nr_bits;
- int q, ind, shift;
- UINT64 C1;
- UINT64 Cstar; // C* represents up to 16 decimal digits ~ 54 bits
- UINT128 fstar;
- UINT128 P128;
- // check for NaN or Infinity
- if ((x & MASK_NAN) == MASK_NAN || (x & MASK_INF) == MASK_INF) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // unpack x
- x_sign = x & MASK_SIGN; // 0 for positive, MASK_SIGN for negative
- // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] =>
- if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) {
- x_exp = (x & MASK_BINARY_EXPONENT2) >> 51; // biased
- C1 = (x & MASK_BINARY_SIG2) | MASK_BINARY_OR2;
- if (C1 > 9999999999999999ull) { // non-canonical
- x_exp = 0;
- C1 = 0;
- }
- } else {
- x_exp = (x & MASK_BINARY_EXPONENT1) >> 53; // biased
- C1 = x & MASK_BINARY_SIG1;
- }
- // check for zeros (possibly from non-canonical values)
- if (C1 == 0x0ull) {
- // x is 0
- res = 0x00000000;
- BID_RETURN (res);
- }
- // x is not special and is not zero
- // q = nr. of decimal digits in x (1 <= q <= 54)
- // determine first the nr. of bits in x
- if (C1 >= 0x0020000000000000ull) { // x >= 2^53
- // split the 64-bit value in two 32-bit halves to avoid rounding errors
- if (C1 >= 0x0000000100000000ull) { // x >= 2^32
- tmp1.d = (double) (C1 >> 32); // exact conversion
- x_nr_bits =
- 33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- } else { // x < 2^32
- tmp1.d = (double) C1; // exact conversion
- x_nr_bits =
- 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- }
- } else { // if x < 2^53
- tmp1.d = (double) C1; // exact conversion
- x_nr_bits =
- 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
- }
- q = nr_digits[x_nr_bits - 1].digits;
- if (q == 0) {
- q = nr_digits[x_nr_bits - 1].digits1;
- if (C1 >= nr_digits[x_nr_bits - 1].threshold_lo)
- q++;
- }
- exp = x_exp - 398; // unbiased exponent
- if ((q + exp) > 10) { // x >= 10^10 ~= 2^33.2... (cannot fit in 32 bits)
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- } else if ((q + exp) == 10) { // x = c(0)c(1)...c(9).c(10)...c(q-1)
- // in this case 2^29.89... ~= 10^9 <= x < 10^10 ~= 2^33.2...
- // so x rounded to an integer may or may not fit in a signed 32-bit int
- // the cases that do not fit are identified here; the ones that fit
- // fall through and will be handled with other cases further,
- // under '1 <= q + exp <= 10'
- if (x_sign) { // if n < 0 and q + exp = 10
- // if n <= -2^31 - 1/2 then n is too large
- // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31+1/2
- // <=> 0.c(0)c(1)...c(q-1) * 10^11 >= 0x500000005, 1<=q<=16
- // <=> C * 10^(11-q) >= 0x500000005, 1<=q<=16
- if (q <= 11) {
- // Note: C * 10^(11-q) has 10 or 11 digits; 0x500000005 has 11 digits
- tmp64 = C1 * ten2k64[11 - q]; // C scaled up to 11-digit int
- // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
- if (tmp64 >= 0x500000005ull) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- } else { // if (q > 11), i.e. 12 <= q <= 16 and so -15 <= exp <= -2
- // C * 10^(11-q) >= 0x500000005 <=>
- // C >= 0x500000005 * 10^(q-11) where 1 <= q - 11 <= 5
- // (scale 2^31+1/2 up)
- // Note: 0x500000005*10^(q-11) has q-1 or q digits, where q <= 16
- tmp64 = 0x500000005ull * ten2k64[q - 11];
- if (C1 >= tmp64) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- }
- } else { // if n > 0 and q + exp = 10
- // if n >= 2^31 - 1/2 then n is too large
- // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31-1/2
- // <=> 0.c(0)c(1)...c(q-1) * 10^11 >= 0x4fffffffb, 1<=q<=16
- // <=> C * 10^(11-q) >= 0x4fffffffb, 1<=q<=16
- if (q <= 11) {
- // Note: C * 10^(11-q) has 10 or 11 digits; 0x500000005 has 11 digits
- tmp64 = C1 * ten2k64[11 - q]; // C scaled up to 11-digit int
- // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
- if (tmp64 >= 0x4fffffffbull) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- } else { // if (q > 11), i.e. 12 <= q <= 16 and so -15 <= exp <= -2
- // C * 10^(11-q) >= 0x4fffffffb <=>
- // C >= 0x4fffffffb * 10^(q-11) where 1 <= q - 11 <= 5
- // (scale 2^31-1/2 up)
- // Note: 0x4fffffffb*10^(q-11) has q-1 or q digits, where q <= 16
- tmp64 = 0x4fffffffbull * ten2k64[q - 11];
- if (C1 >= tmp64) {
- // set invalid flag
- *pfpsf |= INVALID_EXCEPTION;
- // return Integer Indefinite
- res = 0x80000000;
- BID_RETURN (res);
- }
- // else cases that can be rounded to a 32-bit int fall through
- // to '1 <= q + exp <= 10'
- }
- }
- }
- // n is not too large to be converted to int32: -2^31 - 1/2 < n < 2^31 - 1/2
- // Note: some of the cases tested for above fall through to this point
- if ((q + exp) < 0) { // n = +/-0.0...c(0)c(1)...c(q-1)
- // set inexact flag
- *pfpsf |= INEXACT_EXCEPTION;
- // return 0
- res = 0x00000000;
- BID_RETURN (res);
- } else if ((q + exp) == 0) { // n = +/-0.c(0)c(1)...c(q-1)
- // if 0.c(0)c(1)...c(q-1) <= 0.5 <=> c(0)c(1)...c(q-1) <= 5 * 10^(q-1)
- // res = 0
- // else
- // res = +/-1
- ind = q - 1;
- if (C1 < midpoint64[ind]) {
- res = 0x00000000; // return 0
- } else if (x_sign) { // n < 0
- res = 0xffffffff; // return -1
- } else { // n > 0
- res = 0x00000001; // return +1
- }
- // set inexact flag
- *pfpsf |= INEXACT_EXCEPTION;
- } else { // if (1 <= q + exp <= 10, 1 <= q <= 16, -15 <= exp <= 9)
- // -2^31-1/2 <= x <= -1 or 1 <= x < 2^31-1/2 so x can be rounded
- // to nearest away to a 32-bit signed integer
- if (exp < 0) { // 2 <= q <= 16, -15 <= exp <= -1, 1 <= q + exp <= 10
- ind = -exp; // 1 <= ind <= 15; ind is a synonym for 'x'
- // chop off ind digits from the lower part of C1
- // C1 = C1 + 1/2 * 10^ind where the result C1 fits in 64 bits
- C1 = C1 + midpoint64[ind - 1];
- // calculate C* and f*
- // C* is actually floor(C*) in this case
- // C* and f* need shifting and masking, as shown by
- // shiftright128[] and maskhigh128[]
- // 1 <= x <= 15
- // kx = 10^(-x) = ten2mk64[ind - 1]
- // C* = (C1 + 1/2 * 10^x) * 10^(-x)
- // the approximation of 10^(-x) was rounded up to 54 bits
- __mul_64x64_to_128MACH (P128, C1, ten2mk64[ind - 1]);
- Cstar = P128.w[1];
- fstar.w[1] = P128.w[1] & maskhigh128[ind - 1];
- fstar.w[0] = P128.w[0];
- // the top Ex bits of 10^(-x) are T* = ten2mk128trunc[ind].w[0], e.g.
- // if x=1, T*=ten2mk128trunc[0].w[0]=0x1999999999999999
- // C* = floor(C*)-1 (logical right shift; C* has p decimal digits,
- // correct by Pr. 1)
- // n = C* * 10^(e+x)
- // shift right C* by Ex-64 = shiftright128[ind]
- shift = shiftright128[ind - 1]; // 0 <= shift <= 39
- Cstar = Cstar >> shift;
- // determine inexactness of the rounding of C*
- // if (0 < f* - 1/2 < 10^(-x)) then
- // the result is exact
- // else // if (f* - 1/2 > T*) then
- // the result is inexact
- if (ind - 1 <= 2) {
- if (fstar.w[0] > 0x8000000000000000ull) {
- // f* > 1/2 and the result may be exact
- tmp64 = fstar.w[0] - 0x8000000000000000ull; // f* - 1/2
- if ((tmp64 > ten2mk128trunc[ind - 1].w[1])) {
- // ten2mk128trunc[ind -1].w[1] is identical to
- // ten2mk128[ind -1].w[1]
- // set the inexact flag
- *pfpsf |= INEXACT_EXCEPTION;
- } // else the result is exact
- } else { // the result is inexact; f2* <= 1/2
- // set the inexact flag
- *pfpsf |= INEXACT_EXCEPTION;
- }
- } else { // if 3 <= ind - 1 <= 14
- if (fstar.w[1] > onehalf128[ind - 1] ||
- (fstar.w[1] == onehalf128[ind - 1] && fstar.w[0])) {
- // f2* > 1/2 and the result may be exact
- // Calculate f2* - 1/2
- tmp64 = fstar.w[1] - onehalf128[ind - 1];
- if (tmp64 || fstar.w[0] > ten2mk128trunc[ind - 1].w[1]) {
- // ten2mk128trunc[ind -1].w[1] is identical to
- // ten2mk128[ind -1].w[1]
- // set the inexact flag
- *pfpsf |= INEXACT_EXCEPTION;
- } // else the result is exact
- } else { // the result is inexact; f2* <= 1/2
- // set the inexact flag
- *pfpsf |= INEXACT_EXCEPTION;
- }
- }
- // if the result was a midpoint it was rounded away from zero
- if (x_sign)
- res = -Cstar;
- else
- res = Cstar;
- } else if (exp == 0) {
- // 1 <= q <= 10
- // res = +/-C (exact)
- if (x_sign)
- res = -C1;
- else
- res = C1;
- } else { // if (exp > 0) => 1 <= exp <= 9, 1 <= q < 9, 2 <= q + exp <= 10
- // res = +/-C * 10^exp (exact)
- if (x_sign)
- res = -C1 * ten2k64[exp];
- else
- res = C1 * ten2k64[exp];
- }
- }
- BID_RETURN (res);
- }
|