Mercurial > hg > octave-shane > gnulib-hg
annotate lib/exp2.c @ 17282:a4996fb12971
locale: port to Solaris 2.6 and 7 + GNU gettext
* lib/locale.in.h: Just include_next <locale.h> when
being invoked recursively. This prevents problems on Solaris 2.6 and 7
when combining the localename module with GNU gettext 0.18.2.
Problem reported by Tom G. Christensen in
<http://lists.gnu.org/archive/html/bug-gnulib/2013-01/msg00084.html>.
| author | Paul Eggert <eggert@cs.ucla.edu> |
|---|---|
| date | Thu, 10 Jan 2013 13:24:13 -0800 |
| parents | e542fd46ad6f |
| children | 344018b6e5d7 |
| rev | line source |
|---|---|
| 16616 | 1 /* Exponential base 2 function. |
|
17249
e542fd46ad6f
maint: update all copyright year number ranges
Eric Blake <eblake@redhat.com>
parents:
16651
diff
changeset
|
2 Copyright (C) 2012-2013 Free Software Foundation, Inc. |
| 16616 | 3 |
| 4 This program is free software: you can redistribute it and/or modify | |
| 5 it under the terms of the GNU General Public License as published by | |
| 6 the Free Software Foundation; either version 3 of the License, or | |
| 7 (at your option) any later version. | |
| 8 | |
| 9 This program is distributed in the hope that it will be useful, | |
| 10 but WITHOUT ANY WARRANTY; without even the implied warranty of | |
| 11 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the | |
| 12 GNU General Public License for more details. | |
| 13 | |
| 14 You should have received a copy of the GNU General Public License | |
| 15 along with this program. If not, see <http://www.gnu.org/licenses/>. */ | |
| 16 | |
| 17 #include <config.h> | |
| 18 | |
| 19 /* Specification. */ | |
| 20 #include <math.h> | |
| 21 | |
| 22 #include <float.h> | |
| 23 | |
| 24 /* Best possible approximation of log(2) as a 'double'. */ | |
| 25 #define LOG2 0.693147180559945309417232121458176568075 | |
| 26 | |
| 27 /* Best possible approximation of 1/log(2) as a 'double'. */ | |
| 28 #define LOG2_INVERSE 1.44269504088896340735992468100189213743 | |
| 29 | |
| 30 /* Best possible approximation of log(2)/256 as a 'double'. */ | |
| 31 #define LOG2_BY_256 0.00270760617406228636491106297444600221904 | |
| 32 | |
| 33 /* Best possible approximation of 256/log(2) as a 'double'. */ | |
| 34 #define LOG2_BY_256_INVERSE 369.329930467574632284140718336484387181 | |
| 35 | |
| 36 double | |
| 37 exp2 (double x) | |
| 38 { | |
| 39 /* exp2(x) = exp(x*log(2)). | |
| 40 If we would compute it like this, there would be rounding errors for | |
| 41 integer or near-integer values of x. To avoid these, we inline the | |
| 42 algorithm for exp(), and the multiplication with log(2) cancels a | |
| 43 division by log(2). */ | |
| 44 | |
| 45 if (isnand (x)) | |
| 46 return x; | |
| 47 | |
| 48 if (x > (double) DBL_MAX_EXP) | |
| 49 /* x > DBL_MAX_EXP | |
| 50 hence exp2(x) > 2^DBL_MAX_EXP, overflows to Infinity. */ | |
| 51 return HUGE_VAL; | |
| 52 | |
| 53 if (x < (double) (DBL_MIN_EXP - 1 - DBL_MANT_DIG)) | |
| 54 /* x < (DBL_MIN_EXP - 1 - DBL_MANT_DIG) | |
| 55 hence exp2(x) < 2^(DBL_MIN_EXP-1-DBL_MANT_DIG), | |
| 56 underflows to zero. */ | |
| 57 return 0.0; | |
| 58 | |
| 59 /* Decompose x into | |
| 60 x = n + m/256 + y/log(2) | |
| 61 where | |
| 62 n is an integer, | |
| 63 m is an integer, -128 <= m <= 128, | |
| 64 y is a number, |y| <= log(2)/512 + epsilon = 0.00135... | |
| 65 Then | |
| 66 exp2(x) = 2^n * exp(m * log(2)/256) * exp(y) | |
| 67 The first factor is an ldexpl() call. | |
| 68 The second factor is a table lookup. | |
| 69 The third factor is computed | |
| 70 - either as sinh(y) + cosh(y) | |
| 71 where sinh(y) is computed through the power series: | |
| 72 sinh(y) = y + y^3/3! + y^5/5! + ... | |
| 73 and cosh(y) is computed as hypot(1, sinh(y)), | |
| 74 - or as exp(2*z) = (1 + tanh(z)) / (1 - tanh(z)) | |
| 75 where z = y/2 | |
| 76 and tanh(z) is computed through its power series: | |
| 77 tanh(z) = z | |
| 78 - 1/3 * z^3 | |
| 79 + 2/15 * z^5 | |
| 80 - 17/315 * z^7 | |
| 81 + 62/2835 * z^9 | |
| 82 - 1382/155925 * z^11 | |
| 83 + 21844/6081075 * z^13 | |
| 84 - 929569/638512875 * z^15 | |
| 85 + ... | |
| 16651 | 86 Since |z| <= log(2)/1024 < 0.0007, the relative contribution of the |
| 87 z^7 term is < 0.0007^6 < 2^-60 <= 2^-DBL_MANT_DIG, therefore we can | |
| 88 truncate the series after the z^5 term. */ | |
| 16616 | 89 |
| 90 { | |
| 91 double nm = round (x * 256.0); /* = 256 * n + m */ | |
| 92 double z = (x * 256.0 - nm) * (LOG2_BY_256 * 0.5); | |
| 93 | |
| 94 /* Coefficients of the power series for tanh(z). */ | |
| 95 #define TANH_COEFF_1 1.0 | |
| 96 #define TANH_COEFF_3 -0.333333333333333333333333333333333333334 | |
| 97 #define TANH_COEFF_5 0.133333333333333333333333333333333333334 | |
| 98 #define TANH_COEFF_7 -0.053968253968253968253968253968253968254 | |
| 99 #define TANH_COEFF_9 0.0218694885361552028218694885361552028218 | |
| 100 #define TANH_COEFF_11 -0.00886323552990219656886323552990219656886 | |
| 101 #define TANH_COEFF_13 0.00359212803657248101692546136990581435026 | |
| 102 #define TANH_COEFF_15 -0.00145583438705131826824948518070211191904 | |
| 103 | |
| 104 double z2 = z * z; | |
| 105 double tanh_z = | |
| 106 ((TANH_COEFF_5 | |
| 107 * z2 + TANH_COEFF_3) | |
| 108 * z2 + TANH_COEFF_1) | |
| 109 * z; | |
| 110 | |
| 111 double exp_y = (1.0 + tanh_z) / (1.0 - tanh_z); | |
| 112 | |
| 113 int n = (int) round (nm * (1.0 / 256.0)); | |
| 114 int m = (int) nm - 256 * n; | |
| 115 | |
| 116 /* exp_table[i] = exp((i - 128) * log(2)/256). | |
| 117 Computed in GNU clisp through | |
| 118 (setf (long-float-digits) 128) | |
| 119 (setq a 0L0) | |
| 120 (setf (long-float-digits) 256) | |
| 121 (dotimes (i 257) | |
| 122 (format t " ~D,~%" | |
| 123 (float (exp (* (/ (- i 128) 256) (log 2L0))) a))) */ | |
| 124 static const double exp_table[257] = | |
| 125 { | |
| 126 0.707106781186547524400844362104849039284, | |
| 127 0.709023942160207598920563322257676190836, | |
| 128 0.710946301084582779904674297352120049962, | |
| 129 0.71287387205274715340350157671438300618, | |
| 130 0.714806669195985005617532889137569953044, | |
| 131 0.71674470668389442125974978427737336719, | |
| 132 0.71868799872449116280161304224785251353, | |
| 133 0.720636559564312831364255957304947586072, | |
| 134 0.72259040348852331001850312073583545284, | |
| 135 0.724549544821017490259402705487111270714, | |
| 136 0.726513997924526282423036245842287293786, | |
| 137 0.728483777200721910815451524818606761737, | |
| 138 0.730458897090323494325651445155310766577, | |
| 139 0.732439372073202913296664682112279175616, | |
| 140 0.734425216668490963430822513132890712652, | |
| 141 0.736416445434683797507470506133110286942, | |
| 142 0.738413072969749655693453740187024961962, | |
| 143 0.740415113911235885228829945155951253966, | |
| 144 0.742422582936376250272386395864403155277, | |
| 145 0.744435494762198532693663597314273242753, | |
| 146 0.746453864145632424600321765743336770838, | |
| 147 0.748477705883617713391824861712720862423, | |
| 148 0.750507034813212760132561481529764324813, | |
| 149 0.752541865811703272039672277899716132493, | |
| 150 0.75458221379671136988300977551659676571, | |
| 151 0.756628093726304951096818488157633113612, | |
| 152 0.75867952059910734940489114658718937343, | |
| 153 0.760736509454407291763130627098242426467, | |
| 154 0.762799075372269153425626844758470477304, | |
| 155 0.76486723347364351194254345936342587308, | |
| 156 0.766940998920478000900300751753859329456, | |
| 157 0.769020386915828464216738479594307884331, | |
| 158 0.771105412703970411806145931045367420652, | |
| 159 0.773196091570510777431255778146135325272, | |
| 160 0.77529243884249997956151370535341912283, | |
| 161 0.777394469888544286059157168801667390437, | |
| 162 0.779502200118918483516864044737428940745, | |
| 163 0.781615644985678852072965367573877941354, | |
| 164 0.783734819982776446532455855478222575498, | |
| 165 0.78585974064617068462428149076570281356, | |
| 166 0.787990422553943243227635080090952504452, | |
| 167 0.790126881326412263402248482007960521995, | |
| 168 0.79226913262624686505993407346567890838, | |
| 169 0.794417192158581972116898048814333564685, | |
| 170 0.796571075671133448968624321559534367934, | |
| 171 0.798730798954313549131410147104316569576, | |
| 172 0.800896377841346676896923120795476813684, | |
| 173 0.803067828208385462848443946517563571584, | |
| 174 0.805245165974627154089760333678700291728, | |
| 175 0.807428407102430320039984581575729114268, | |
| 176 0.809617567597431874649880866726368203972, | |
| 177 0.81181266350866441589760797777344082227, | |
| 178 0.814013710928673883424109261007007338614, | |
| 179 0.816220725993637535170713864466769240053, | |
| 180 0.818433724883482243883852017078007231025, | |
| 181 0.82065272382200311435413206848451310067, | |
| 182 0.822877739076982422259378362362911222833, | |
| 183 0.825108786960308875483586738272485101678, | |
| 184 0.827345883828097198786118571797909120834, | |
| 185 0.829589046080808042697824787210781231927, | |
| 186 0.831838290163368217523168228488195222638, | |
| 187 0.834093632565291253329796170708536192903, | |
| 188 0.836355089820798286809404612069230711295, | |
| 189 0.83862267850893927589613232455870870518, | |
| 190 0.84089641525371454303112547623321489504, | |
| 191 0.84317631672419664796432298771385230143, | |
| 192 0.84546239963465259098692866759361830709, | |
| 193 0.84775468074466634749045860363936420312, | |
| 194 0.850053176859261734750681286748751167545, | |
| 195 0.852357904829025611837203530384718316326, | |
| 196 0.854668881550231413551897437515331498025, | |
| 197 0.856986123964963019301812477839166009452, | |
| 198 0.859309649061238957814672188228156252257, | |
| 199 0.861639473873136948607517116872358729753, | |
| 200 0.863975615480918781121524414614366207052, | |
| 201 0.866318091011155532438509953514163469652, | |
| 202 0.868666917636853124497101040936083380124, | |
| 203 0.871022112577578221729056715595464682243, | |
| 204 0.873383693099584470038708278290226842228, | |
| 205 0.875751676515939078050995142767930296012, | |
| 206 0.878126080186649741556080309687656610647, | |
| 207 0.880506921518791912081045787323636256171, | |
| 208 0.882894217966636410521691124969260937028, | |
| 209 0.885287987031777386769987907431242017412, | |
| 210 0.88768824626326062627527960009966160388, | |
| 211 0.89009501325771220447985955243623523504, | |
| 212 0.892508305659467490072110281986409916153, | |
| 213 0.8949281411607004980029443898876582985, | |
| 214 0.897354537501553593213851621063890907178, | |
| 215 0.899787512470267546027427696662514569756, | |
| 216 0.902227083903311940153838631655504844215, | |
| 217 0.904673269685515934269259325789226871994, | |
| 218 0.907126087750199378124917300181170171233, | |
| 219 0.909585556079304284147971563828178746372, | |
| 220 0.91205169270352665549806275316460097744, | |
| 221 0.914524515702448671545983912696158354092, | |
| 222 0.91700404320467123174354159479414442804, | |
| 223 0.919490293387946858856304371174663918816, | |
| 224 0.921983284479312962533570386670938449637, | |
| 225 0.92448303475522546419252726694739603678, | |
| 226 0.92698956254169278419622653516884831976, | |
| 227 0.929502886214410192307650717745572682403, | |
| 228 0.932023024198894522404814545597236289343, | |
| 229 0.934549994970619252444512104439799143264, | |
| 230 0.93708381705514995066499947497722326722, | |
| 231 0.93962450902828008902058735120448448827, | |
| 232 0.942172089516167224843810351983745154882, | |
| 233 0.944726577195469551733539267378681531548, | |
| 234 0.947287990793482820670109326713462307376, | |
| 235 0.949856349088277632361251759806996099924, | |
| 236 0.952431670908837101825337466217860725517, | |
| 237 0.955013975135194896221170529572799135168, | |
| 238 0.957603280698573646936305635147915443924, | |
| 239 0.960199606581523736948607188887070611744, | |
| 240 0.962802971818062464478519115091191368377, | |
| 241 0.965413395493813583952272948264534783197, | |
| 242 0.968030896746147225299027952283345762418, | |
| 243 0.970655494764320192607710617437589705184, | |
| 244 0.973287208789616643172102023321302921373, | |
| 245 0.97592605811548914795551023340047499377, | |
| 246 0.978572062087700134509161125813435745597, | |
| 247 0.981225240104463713381244885057070325016, | |
| 248 0.983885611616587889056366801238014683926, | |
| 249 0.98655319612761715646797006813220671315, | |
| 250 0.989228013193975484129124959065583667775, | |
| 251 0.99191008242510968492991311132615581644, | |
| 252 0.994599423483633175652477686222166314457, | |
| 253 0.997296056085470126257659913847922601123, | |
| 254 1.0, | |
| 255 1.00271127505020248543074558845036204047, | |
| 256 1.0054299011128028213513839559347998147, | |
| 257 1.008155898118417515783094890817201039276, | |
| 258 1.01088928605170046002040979056186052439, | |
| 259 1.013630084951489438840258929063939929597, | |
| 260 1.01637831491095303794049311378629406276, | |
| 261 1.0191339960777379496848780958207928794, | |
| 262 1.02189714865411667823448013478329943978, | |
| 263 1.02466779289713564514828907627081492763, | |
| 264 1.0274459491187636965388611939222137815, | |
| 265 1.030231637686041012871707902453904567093, | |
| 266 1.033024879021228422500108283970460918086, | |
| 267 1.035825693601957120029983209018081371844, | |
| 268 1.03863410196137879061243669795463973258, | |
| 269 1.04145012468831614126454607901189312648, | |
| 270 1.044273782427413840321966478739929008784, | |
| 271 1.04710509587928986612990725022711224056, | |
| 272 1.04994408580068726608203812651590790906, | |
| 273 1.05279077300462632711989120298074630319, | |
| 274 1.05564517836055715880834132515293865216, | |
| 275 1.058507322794512690105772109683716645074, | |
| 276 1.061377227289262080950567678003883726294, | |
| 277 1.06425491288446454978861125700158022068, | |
| 278 1.06714040067682361816952112099280916261, | |
| 279 1.0700337118202417735424119367576235685, | |
| 280 1.072934867525975551385035450873827585343, | |
| 281 1.075843889062791037803228648476057074063, | |
| 282 1.07876079775711979374068003743848295849, | |
| 283 1.081685614993215201942115594422531125643, | |
| 284 1.08461836221330923781610517190661434161, | |
| 285 1.087559060917769665346797830944039707867, | |
| 286 1.09050773266525765920701065576070797899, | |
| 287 1.09346439907288585422822014625044716208, | |
| 288 1.096429081816376823386138295859248481766, | |
| 289 1.09940180263022198546369696823882990404, | |
| 290 1.10238258330784094355641420942564685751, | |
| 291 1.10537144570174125558827469625695031104, | |
| 292 1.108368411723678638009423649426619850137, | |
| 293 1.111373503344817603850149254228916637444, | |
| 294 1.1143867425958925363088129569196030678, | |
| 295 1.11740815156736919905457996308578026665, | |
| 296 1.12043775240960668442900387986631301277, | |
| 297 1.123475567333019800733729739775321431954, | |
| 298 1.12652161860824189979479864378703477763, | |
| 299 1.129575928566288145997264988840249825907, | |
| 300 1.13263851959871922798707372367762308438, | |
| 301 1.13570941415780551424039033067611701343, | |
| 302 1.13878863475669165370383028384151125472, | |
| 303 1.14187620396956162271229760828788093894, | |
| 304 1.14497214443180421939441388822291589579, | |
| 305 1.14807647884017900677879966269734268003, | |
| 306 1.15118922995298270581775963520198253612, | |
| 307 1.154310420590216039548221528724806960684, | |
| 308 1.157440073633751029613085766293796821106, | |
| 309 1.16057821202749874636945947257609098625, | |
| 310 1.16372485877757751381357359909218531234, | |
| 311 1.166880036952481570555516298414089287834, | |
| 312 1.170043769683250188080259035792738573, | |
| 313 1.17321608016363724753480435451324538889, | |
| 314 1.176396991650281276284645728483848641054, | |
| 315 1.17958652746287594548610056676944051898, | |
| 316 1.182784710984341029924457204693850757966, | |
| 317 1.18599156566099383137126564953421556374, | |
| 318 1.18920711500272106671749997056047591529, | |
| 319 1.19243138258315122214272755814543101148, | |
| 320 1.195664392039827374583837049865451975705, | |
| 321 1.19890616707438048177030255797630020695, | |
| 322 1.202156731452703142096396957497765876003, | |
| 323 1.205416109005123825604211432558411335666, | |
| 324 1.208684323626581577354792255889216998484, | |
| 325 1.21196139927680119446816891773249304545, | |
| 326 1.215247359980468878116520251338798457624, | |
| 327 1.218542229827408361758207148117394510724, | |
| 328 1.221846032972757516903891841911570785836, | |
| 329 1.225158793637145437709464594384845353707, | |
| 330 1.22848053610687000569400895779278184036, | |
| 331 1.2318112847340759358845566532127948166, | |
| 332 1.235151063936933305692912507415415760294, | |
| 333 1.238499898199816567833368865859612431545, | |
| 334 1.24185781207348404859367746872659560551, | |
| 335 1.24522483017525793277520496748615267417, | |
| 336 1.24860097718920473662176609730249554519, | |
| 337 1.25198627786631627006020603178920359732, | |
| 338 1.255380757024691089579390657442301194595, | |
| 339 1.25878443954971644307786044181516261876, | |
| 340 1.26219735039425070801401025851841645967, | |
| 341 1.265619514578806324196273999873453036296, | |
| 342 1.26905095719173322255441908103233800472, | |
| 343 1.27249170338940275123669204418460217677, | |
| 344 1.27594177839639210038120243475928938891, | |
| 345 1.27940120750566922691358797002785254596, | |
| 346 1.28287001607877828072666978102151405111, | |
| 347 1.286348229546025533601482208069738348355, | |
| 348 1.28983587340666581223274729549155218968, | |
| 349 1.293332973229089436725559789048704304684, | |
| 350 1.296839554651009665933754117792451159835, | |
| 351 1.30035564337965065101414056707091779129, | |
| 352 1.30388126519193589857452364895199736833, | |
| 353 1.30741644593467724479715157747196172848, | |
| 354 1.310961211524764341922991786330755849366, | |
| 355 1.314515587949354658485983613383997794965, | |
| 356 1.318079601266063994690185647066116617664, | |
| 357 1.32165327760315751432651181233060922616, | |
| 358 1.32523664315974129462953709549872167411, | |
| 359 1.32882972420595439547865089632866510792, | |
| 360 1.33243254708316144935164337949073577407, | |
| 361 1.33604513820414577344262790437186975929, | |
| 362 1.33966752405330300536003066972435257602, | |
| 363 1.34329973118683526382421714618163087542, | |
| 364 1.346941786232945835788173713229537282075, | |
| 365 1.35059371589203439140852219606013396004, | |
| 366 1.35425554693689272829801474014070280434, | |
| 367 1.357927306212901046494536695671766697446, | |
| 368 1.36160902063822475558553593883194147464, | |
| 369 1.36530071720401181543069836033754285543, | |
| 370 1.36900242297459061192960113298219283217, | |
| 371 1.37271416508766836928499785714471721579, | |
| 372 1.37643597075453010021632280551868696026, | |
| 373 1.380167867260238095581945274358283464697, | |
| 374 1.383909881963831954872659527265192818, | |
| 375 1.387662042298529159042861017950775988896, | |
| 376 1.39142437577192618714983552956624344668, | |
| 377 1.395196909966200178275574599249220994716, | |
| 378 1.398979672538311140209528136715194969206, | |
| 379 1.40277269122020470637471352433337881711, | |
| 380 1.40657599381901544248361973255451684411, | |
| 381 1.410389608217270704414375128268675481145, | |
| 382 1.41421356237309504880168872420969807857 | |
| 383 }; | |
| 384 | |
| 385 return ldexp (exp_table[128 + m] * exp_y, n); | |
| 386 } | |
| 387 } |