1.file "scalbf.s" 2 3 4// Copyright (c) 2000 - 2003, Intel Corporation 5// All rights reserved. 6// 7// 8// Redistribution and use in source and binary forms, with or without 9// modification, are permitted provided that the following conditions are 10// met: 11// 12// * Redistributions of source code must retain the above copyright 13// notice, this list of conditions and the following disclaimer. 14// 15// * Redistributions in binary form must reproduce the above copyright 16// notice, this list of conditions and the following disclaimer in the 17// documentation and/or other materials provided with the distribution. 18// 19// * The name of Intel Corporation may not be used to endorse or promote 20// products derived from this software without specific prior written 21// permission. 22 23// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS 24// "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT 25// LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR 26// A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL INTEL OR ITS 27// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, 28// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, 29// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR 30// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY 31// OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY OR TORT (INCLUDING 32// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS 33// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 34// 35// Intel Corporation is the author of this code, and requests that all 36// problem reports or change requests be submitted to it directly at 37// http://www.intel.com/software/products/opensource/libraries/num.htm. 38// 39// History 40//============================================================== 41// 02/02/00 Initial version 42// 01/26/01 Scalb completely reworked and now standalone version 43// 05/20/02 Cleaned up namespace and sf0 syntax 44// 02/10/03 Reordered header: .section, .global, .proc, .align 45// 08/06/03 Improved performance 46// 47// API 48//============================================================== 49// float = scalbf (float x, float n) 50// input floating point f8 and floating point f9 51// output floating point f8 52// 53// int_type = 0 if int is 32 bits 54// int_type = 1 if int is 64 bits 55// 56// Returns x* 2**n using an fma and detects overflow 57// and underflow. 58// 59// 60// Strategy: 61// Compute biased exponent of result exp_Result = N + exp_X 62// Break into ranges: 63// exp_Result > 0x1007e -> Certain overflow 64// exp_Result = 0x1007e -> Possible overflow 65// 0x0ff81 <= exp_Result < 0x1007e -> No over/underflow (main path) 66// 0x0ff81 - 23 <= exp_Result < 0x0ff81 -> Possible underflow 67// exp_Result < 0x0ff81 - 23 -> Certain underflow 68 69FR_Big = f6 70FR_NBig = f7 71FR_Floating_X = f8 72FR_Result = f8 73FR_Floating_N = f9 74FR_Result2 = f9 75FR_Result3 = f10 76FR_Norm_X = f11 77FR_Two_N = f12 78FR_N_float_int = f13 79FR_Norm_N = f14 80 81GR_neg_ov_limit= r14 82GR_big_exp = r14 83GR_N_Biased = r15 84GR_Big = r16 85GR_exp_Result = r18 86GR_pos_ov_limit= r19 87GR_exp_sure_ou = r19 88GR_Bias = r20 89GR_N_as_int = r21 90GR_signexp_X = r22 91GR_exp_X = r23 92GR_exp_mask = r24 93GR_max_exp = r25 94GR_min_exp = r26 95GR_min_den_exp = r27 96GR_Scratch = r28 97GR_signexp_N = r29 98GR_exp_N = r30 99 100GR_SAVE_B0 = r32 101GR_SAVE_GP = r33 102GR_SAVE_PFS = r34 103GR_Parameter_X = r35 104GR_Parameter_Y = r36 105GR_Parameter_RESULT = r37 106GR_Tag = r38 107 108.section .text 109GLOBAL_IEEE754_ENTRY(scalbf) 110 111// 112// Is x NAN, INF, ZERO, +-? 113// Build the exponent Bias 114// 115{ .mfi 116 getf.exp GR_signexp_N = FR_Floating_N // Get signexp of n 117 fclass.m p6,p0 = FR_Floating_X, 0xe7 // @snan | @qnan | @inf | @zero 118 mov GR_Bias = 0x0ffff 119} 120{ .mfi 121 mov GR_Big = 35000 // If N this big then certain overflow 122 fcvt.fx.trunc.s1 FR_N_float_int = FR_Floating_N // Get N in significand 123 nop.i 0 124} 125;; 126 127{ .mfi 128 getf.exp GR_signexp_X = FR_Floating_X // Get signexp of x 129 fclass.m p7,p0 = FR_Floating_N, 0x0b // Test for n=unorm 130 nop.i 0 131} 132// 133// Normalize n 134// 135{ .mfi 136 mov GR_exp_mask = 0x1ffff // Exponent mask 137 fnorm.s1 FR_Norm_N = FR_Floating_N 138 nop.i 0 139} 140;; 141 142// 143// Is n NAN, INF, ZERO, +-? 144// 145{ .mfi 146 mov GR_big_exp = 0x1003e // Exponent at which n is integer 147 fclass.m p9,p0 = FR_Floating_N, 0xe7 // @snan | @qnan | @inf | @zero 148 mov GR_max_exp = 0x1007e // Exponent of maximum float 149} 150// 151// Normalize x 152// 153{ .mfb 154 nop.m 0 155 fnorm.s1 FR_Norm_X = FR_Floating_X 156(p7) br.cond.spnt SCALBF_N_UNORM // Branch if n=unorm 157} 158;; 159 160SCALBF_COMMON1: 161// Main path continues. Also return here from u=unorm path. 162// Handle special cases if x = Nan, Inf, Zero 163{ .mfb 164 nop.m 0 165 fcmp.lt.s1 p7,p0 = FR_Floating_N, f0 // Test N negative 166(p6) br.cond.spnt SCALBF_NAN_INF_ZERO 167} 168;; 169 170// Handle special cases if n = Nan, Inf, Zero 171{ .mfi 172 getf.sig GR_N_as_int = FR_N_float_int // Get n from significand 173 fclass.m p8,p0 = FR_Floating_X, 0x0b // Test for x=unorm 174 mov GR_exp_sure_ou = 0x1000e // Exp_N where x*2^N sure over/under 175} 176{ .mfb 177 mov GR_min_exp = 0x0ff81 // Exponent of minimum float 178 fcvt.xf FR_N_float_int = FR_N_float_int // Convert N to FP integer 179(p9) br.cond.spnt SCALBF_NAN_INF_ZERO 180} 181;; 182 183{ .mmi 184 and GR_exp_N = GR_exp_mask, GR_signexp_N // Get exponent of N 185(p7) sub GR_Big = r0, GR_Big // Limit for N 186 nop.i 0 187} 188;; 189 190{ .mib 191 cmp.lt p9,p0 = GR_exp_N, GR_big_exp // N possible non-integer? 192 cmp.ge p6,p0 = GR_exp_N, GR_exp_sure_ou // N certain over/under? 193(p8) br.cond.spnt SCALBF_X_UNORM // Branch if x=unorm 194} 195;; 196 197SCALBF_COMMON2: 198// Main path continues. Also return here from x=unorm path. 199// Create biased exponent for 2**N 200{ .mmi 201(p6) mov GR_N_as_int = GR_Big // Limit N 202;; 203 add GR_N_Biased = GR_Bias,GR_N_as_int 204 nop.i 0 205} 206;; 207 208{ .mfi 209 setf.exp FR_Two_N = GR_N_Biased // Form 2**N 210(p9) fcmp.neq.unc.s1 p9,p0 = FR_Norm_N, FR_N_float_int // Test if N an integer 211 and GR_exp_X = GR_exp_mask, GR_signexp_X // Get exponent of X 212} 213;; 214 215// 216// Compute biased result exponent 217// Branch if N is not an integer 218// 219{ .mib 220 add GR_exp_Result = GR_exp_X, GR_N_as_int 221 mov GR_min_den_exp = 0x0ff81 - 23 // Exponent of min denorm float 222(p9) br.cond.spnt SCALBF_N_NOT_INT 223} 224;; 225 226// 227// Raise Denormal operand flag with compare 228// Do final operation 229// 230{ .mfi 231 cmp.lt p7,p6 = GR_exp_Result, GR_max_exp // Test no overflow 232 fcmp.ge.s0 p0,p11 = FR_Floating_X,FR_Floating_N // Dummy to set denorm 233 cmp.lt p9,p0 = GR_exp_Result, GR_min_den_exp // Test sure underflow 234} 235{ .mfb 236 nop.m 0 237 fma.s.s0 FR_Result = FR_Two_N,FR_Norm_X,f0 238(p9) br.cond.spnt SCALBF_UNDERFLOW // Branch if certain underflow 239} 240;; 241 242{ .mib 243(p6) cmp.gt.unc p6,p8 = GR_exp_Result, GR_max_exp // Test sure overflow 244(p7) cmp.ge.unc p7,p9 = GR_exp_Result, GR_min_exp // Test no over/underflow 245(p7) br.ret.sptk b0 // Return from main path 246} 247;; 248 249{ .bbb 250(p6) br.cond.spnt SCALBF_OVERFLOW // Branch if certain overflow 251(p8) br.cond.spnt SCALBF_POSSIBLE_OVERFLOW // Branch if possible overflow 252(p9) br.cond.spnt SCALBF_POSSIBLE_UNDERFLOW // Branch if possible underflow 253} 254;; 255 256// Here if possible underflow. 257// Resulting exponent: 0x0ff81-23 <= exp_Result < 0x0ff81 258SCALBF_POSSIBLE_UNDERFLOW: 259// 260// Here if possible overflow. 261// Resulting exponent: 0x1007e = exp_Result 262SCALBF_POSSIBLE_OVERFLOW: 263 264// Set up necessary status fields 265// 266// S0 user supplied status 267// S2 user supplied status + WRE + TD (Overflows) 268// S3 user supplied status + FZ + TD (Underflows) 269// 270{ .mfi 271 mov GR_pos_ov_limit = 0x1007f // Exponent for positive overflow 272 fsetc.s3 0x7F,0x41 273 nop.i 0 274} 275{ .mfi 276 mov GR_neg_ov_limit = 0x3007f // Exponent for negative overflow 277 fsetc.s2 0x7F,0x42 278 nop.i 0 279} 280;; 281 282// 283// Do final operation with s2 and s3 284// 285{ .mfi 286 setf.exp FR_NBig = GR_neg_ov_limit 287 fma.s.s3 FR_Result3 = FR_Two_N,FR_Norm_X,f0 288 nop.i 0 289} 290{ .mfi 291 setf.exp FR_Big = GR_pos_ov_limit 292 fma.s.s2 FR_Result2 = FR_Two_N,FR_Norm_X,f0 293 nop.i 0 294} 295;; 296 297// Check for overflow or underflow. 298// Restore s3 299// Restore s2 300// 301{ .mfi 302 nop.m 0 303 fsetc.s3 0x7F,0x40 304 nop.i 0 305} 306{ .mfi 307 nop.m 0 308 fsetc.s2 0x7F,0x40 309 nop.i 0 310} 311;; 312 313// 314// Is the result zero? 315// 316{ .mfi 317 nop.m 0 318 fclass.m p6, p0 = FR_Result3, 0x007 319 nop.i 0 320} 321{ .mfi 322 nop.m 0 323 fcmp.ge.s1 p7, p8 = FR_Result2 , FR_Big 324 nop.i 0 325} 326;; 327 328// 329// Detect masked underflow - Tiny + Inexact Only 330// 331{ .mfi 332 nop.m 0 333(p6) fcmp.neq.unc.s1 p6, p0 = FR_Result , FR_Result2 334 nop.i 0 335} 336;; 337 338// 339// Is result bigger the allowed range? 340// Branch out for underflow 341// 342{ .mfb 343 nop.m 0 344(p8) fcmp.le.unc.s1 p9, p10 = FR_Result2 , FR_NBig 345(p6) br.cond.spnt SCALBF_UNDERFLOW 346} 347;; 348 349// 350// Branch out for overflow 351// 352{ .bbb 353(p7) br.cond.spnt SCALBF_OVERFLOW 354(p9) br.cond.spnt SCALBF_OVERFLOW 355 br.ret.sptk b0 // Return from main path. 356} 357;; 358 359// Here if result overflows 360SCALBF_OVERFLOW: 361{ .mib 362 alloc r32=ar.pfs,3,0,4,0 363 addl GR_Tag = 55, r0 // Set error tag for overflow 364 br.cond.sptk __libm_error_region // Call error support for overflow 365} 366;; 367 368// Here if result underflows 369SCALBF_UNDERFLOW: 370{ .mib 371 alloc r32=ar.pfs,3,0,4,0 372 addl GR_Tag = 56, r0 // Set error tag for underflow 373 br.cond.sptk __libm_error_region // Call error support for underflow 374} 375;; 376 377SCALBF_NAN_INF_ZERO: 378 379// 380// Before entry, N has been converted to a fp integer in significand of 381// FR_N_float_int 382// 383// Convert N_float_int to floating point value 384// 385{ .mfi 386 getf.sig GR_N_as_int = FR_N_float_int 387 fclass.m p6,p0 = FR_Floating_N, 0xc3 //@snan | @qnan 388 nop.i 0 389} 390{ .mfi 391 addl GR_Scratch = 1,r0 392 fcvt.xf FR_N_float_int = FR_N_float_int 393 nop.i 0 394} 395;; 396 397{ .mfi 398 nop.m 0 399 fclass.m p7,p0 = FR_Floating_X, 0xc3 //@snan | @qnan 400 shl GR_Scratch = GR_Scratch,63 401} 402;; 403 404{ .mfi 405 nop.m 0 406 fclass.m p8,p0 = FR_Floating_N, 0x21 // @inf 407 nop.i 0 408} 409{ .mfi 410 nop.m 0 411 fclass.m p9,p0 = FR_Floating_N, 0x22 // @-inf 412 nop.i 0 413} 414;; 415 416// 417// Either X or N is a Nan, return result and possible raise invalid. 418// 419{ .mfb 420 nop.m 0 421(p6) fma.s.s0 FR_Result = FR_Floating_N,FR_Floating_X,f0 422(p6) br.ret.spnt b0 423} 424;; 425 426{ .mfb 427 nop.m 0 428(p7) fma.s.s0 FR_Result = FR_Floating_N,FR_Floating_X,f0 429(p7) br.ret.spnt b0 430} 431;; 432 433// 434// If N + Inf do something special 435// For N = -Inf, create Int 436// 437{ .mfb 438 nop.m 0 439(p8) fma.s.s0 FR_Result = FR_Floating_X, FR_Floating_N,f0 440(p8) br.ret.spnt b0 441} 442{ .mfi 443 nop.m 0 444(p9) fnma.s.s0 FR_Floating_N = FR_Floating_N, f1, f0 445 nop.i 0 446} 447;; 448 449// 450// If N==-Inf,return x/(-N) 451// 452{ .mfb 453 cmp.ne p7,p0 = GR_N_as_int,GR_Scratch 454(p9) frcpa.s0 FR_Result,p0 = FR_Floating_X,FR_Floating_N 455(p9) br.ret.spnt b0 456} 457;; 458 459// 460// Is N an integer. 461// 462{ .mfi 463 nop.m 0 464(p7) fcmp.neq.unc.s1 p7,p0 = FR_Norm_N, FR_N_float_int 465 nop.i 0 466} 467;; 468 469// 470// If N not an int, return NaN and raise invalid. 471// 472{ .mfb 473 nop.m 0 474(p7) frcpa.s0 FR_Result,p0 = f0,f0 475(p7) br.ret.spnt b0 476} 477;; 478 479// 480// Always return x in other path. 481// 482{ .mfb 483 nop.m 0 484 fma.s.s0 FR_Result = FR_Floating_X,f1,f0 485 br.ret.sptk b0 486} 487;; 488 489// Here if n not int 490// Return NaN and raise invalid. 491SCALBF_N_NOT_INT: 492{ .mfb 493 nop.m 0 494 frcpa.s0 FR_Result,p0 = f0,f0 495 br.ret.sptk b0 496} 497;; 498 499// Here if n=unorm 500SCALBF_N_UNORM: 501{ .mfb 502 getf.exp GR_signexp_N = FR_Norm_N // Get signexp of normalized n 503 fcvt.fx.trunc.s1 FR_N_float_int = FR_Norm_N // Get N in significand 504 br.cond.sptk SCALBF_COMMON1 // Return to main path 505} 506;; 507 508// Here if x=unorm 509SCALBF_X_UNORM: 510{ .mib 511 getf.exp GR_signexp_X = FR_Norm_X // Get signexp of normalized x 512 nop.i 0 513 br.cond.sptk SCALBF_COMMON2 // Return to main path 514} 515;; 516 517GLOBAL_IEEE754_END(scalbf) 518LOCAL_LIBM_ENTRY(__libm_error_region) 519 520// 521// Get stack address of N 522// 523.prologue 524{ .mfi 525 add GR_Parameter_Y=-32,sp 526 nop.f 0 527.save ar.pfs,GR_SAVE_PFS 528 mov GR_SAVE_PFS=ar.pfs 529} 530// 531// Adjust sp 532// 533{ .mfi 534.fframe 64 535 add sp=-64,sp 536 nop.f 0 537 mov GR_SAVE_GP=gp 538};; 539 540// 541// Store N on stack in correct position 542// Locate the address of x on stack 543// 544{ .mmi 545 stfs [GR_Parameter_Y] = FR_Norm_N,16 546 add GR_Parameter_X = 16,sp 547.save b0, GR_SAVE_B0 548 mov GR_SAVE_B0=b0 549};; 550 551// 552// Store x on the stack. 553// Get address for result on stack. 554// 555.body 556{ .mib 557 stfs [GR_Parameter_X] = FR_Norm_X 558 add GR_Parameter_RESULT = 0,GR_Parameter_Y 559 nop.b 0 560} 561{ .mib 562 stfs [GR_Parameter_Y] = FR_Result 563 add GR_Parameter_Y = -16,GR_Parameter_Y 564 br.call.sptk b0=__libm_error_support# 565};; 566 567// 568// Get location of result on stack 569// 570{ .mmi 571 add GR_Parameter_RESULT = 48,sp 572 nop.m 0 573 nop.i 0 574};; 575 576// 577// Get the new result 578// 579{ .mmi 580 ldfs FR_Result = [GR_Parameter_RESULT] 581.restore sp 582 add sp = 64,sp 583 mov b0 = GR_SAVE_B0 584};; 585 586// 587// Restore gp, ar.pfs and return 588// 589{ .mib 590 mov gp = GR_SAVE_GP 591 mov ar.pfs = GR_SAVE_PFS 592 br.ret.sptk b0 593};; 594 595LOCAL_LIBM_END(__libm_error_region) 596 597.type __libm_error_support#,@function 598.global __libm_error_support# 599