sysdeps/ia64/fpu/e_sqrtf.S

   1 .file "sqrtf.s"
   2
   3
   4 // Copyright (c) 2000 - 2003, Intel Corporation
   5 // All rights reserved.
   6 //
   7 // Contributed 2000 by the Intel Numerics Group, Intel Corporation
   8 //
   9 // Redistribution and use in source and binary forms, with or without
  10 // modification, are permitted provided that the following conditions are
  11 // met:
  12 //
  13 // * Redistributions of source code must retain the above copyright
  14 // notice, this list of conditions and the following disclaimer.
  15 //
  16 // * Redistributions in binary form must reproduce the above copyright
  17 // notice, this list of conditions and the following disclaimer in the
  18 // documentation and/or other materials provided with the distribution.
  19 //
  20 // * The name of Intel Corporation may not be used to endorse or promote
  21 // products derived from this software without specific prior written
  22 // permission.
  23
  24 // THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
  25 // "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
  26 // LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
  27 // A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL INTEL OR ITS
  28 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
  29 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
  30 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
  31 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY
  32 // OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY OR TORT (INCLUDING
  33 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
  34 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
  35 //
  36 // Intel Corporation is the author of this code, and requests that all
  37 // problem reports or change requests be submitted to it directly at
  38 // http://www.intel.com/software/products/opensource/libraries/num.htm.
  39 //
  40 //*********************************************************************
  41 // History:
  42 //
  43 // 02/02/00 Initial version
  44 // 04/04/00 Unwind support added
  45 // 08/15/00 Bundle added after call to __libm_error_support to properly
  46 //          set [the previously overwritten] GR_Parameter_RESULT.
  47 // 05/20/02 Cleaned up namespace and sf0 syntax
  48 // 02/10/03 Reordered header: .section, .global, .proc, .align
  49 //
  50 //*********************************************************************
  51 //
  52 // Function:   Combined sqrtf(x), where
  53 //                         _
  54 //             sqrtf(x) = |x, for single precision x values
  55 //
  56 //********************************************************************
  57 //
  58 // Accuracy:       Correctly Rounded
  59 //
  60 //********************************************************************
  61 //
  62 // Resources Used:
  63 //
  64 //    Floating-Point Registers: f8  (Input and Return Value)
  65 //                              f7 -f14
  66 //
  67 //    General Purpose Registers:
  68 //      r32-r36 (Locals)
  69 //      r37-r40 (Used to pass arguments to error handling routine)
  70 //
  71 //    Predicate Registers:      p6, p7, p8
  72 //
  73 //********************************************************************
  74 //
  75 // IEEE Special Conditions:
  76 //
  77 //    All faults and exceptions should be raised correctly.
  78 //    sqrtf(QNaN) = QNaN
  79 //    sqrtf(SNaN) = QNaN
  80 //    sqrtf(+/-0) = +/-0
  81 //    sqrtf(negative) = QNaN and error handling is called
  82 //
  83 //********************************************************************
  84 //
  85 // Implementation:
  86 //
  87 //  Modified Newton-Raphson Algorithm
  88 //
  89 //********************************************************************
  90
  91
  92 GR_SAVE_B0                    = r34
  93 GR_SAVE_PFS                   = r33
  94 GR_SAVE_GP                    = r35
  95
  96 GR_Parameter_X                = r37
  97 GR_Parameter_Y                = r38
  98 GR_Parameter_RESULT           = r39
  99 GR_Parameter_TAG              = r40
 100
 101 FR_X             = f13
 102 FR_Y             = f0
 103 FR_RESULT        = f8
 104
 105
 106 .section .text
 107 GLOBAL_IEEE754_ENTRY(sqrtf)
 108 { .mlx
 109   // BEGIN SINGLE PRECISION MINIMUM LATENCY SQUARE ROOT ALGORITHM
 110   alloc r32= ar.pfs,0,5,4,0
 111   // exponent of +1/2 in r2
 112   movl r2 = 0x0fffe
 113 } { .mfi
 114   // +1/2 in f12
 115   nop.m 0
 116   frsqrta.s0 f7,p6=f8
 117   nop.i 0;;
 118 } { .mfi
 119   setf.exp f12 = r2
 120   // Step (1)
 121   // y0 = 1/sqrt(a) in f7
 122   fclass.m.unc p7,p8 = f8,0x3A
 123   nop.i 0
 124 } { .mfi
 125   nop.m 0
 126   // Make a copy of x just in case
 127   mov f13 = f8
 128   nop.i 0;;
 129 } { .mfi
 130   nop.m 0
 131   // Step (2)
 132   // H0 = 1/2 * y0 in f9
 133   (p6) fma.s1 f9=f12,f7,f0
 134   nop.i 0
 135 } { .mfi
 136   nop.m 0
 137   // Step (3)
 138   // S0 = a * y0 in f7
 139   (p6) fma.s1 f7=f8,f7,f0
 140   nop.i 0;;
 141 } { .mfi
 142   nop.m 0
 143   // Step (4)
 144   // d = 1/2 - S0 * H0 in f10
 145   (p6) fnma.s1 f10=f7,f9,f12
 146   nop.i 0
 147 } { .mfi
 148   nop.m 0
 149   // Step (0'')
 150   // 3/2 = 1 + 1/2 in f12
 151   (p6) fma.s1 f12=f12,f1,f1
 152   nop.i 0;;
 153 } { .mfi
 154   nop.m 0
 155   // Step (5)
 156   // e = 1 + 3/2 * d in f12
 157   (p6) fma.s1 f12=f12,f10,f1
 158   nop.i 0
 159 } { .mfi
 160   nop.m 0
 161   // Step (6)
 162   // T0 = d * S0 in f11
 163   (p6) fma.s1 f11=f10,f7,f0
 164   nop.i 0;;
 165 } { .mfi
 166   nop.m 0
 167   // Step (7)
 168   // G0 = d * H0 in f10
 169   (p6) fma.s1 f10=f10,f9,f0
 170   nop.i 0;;
 171 } { .mfi
 172   nop.m 0
 173   // Step (8)
 174   // S1 = S0 + e * T0 in f7
 175   (p6) fma.s.s1 f7=f12,f11,f7
 176   nop.i 0;;
 177 } { .mfi
 178   nop.m 0
 179   // Step (9)
 180   // H1 = H0 + e * G0 in f12
 181   (p6) fma.s1 f12=f12,f10,f9
 182   nop.i 0;;
 183 } { .mfi
 184   nop.m 0
 185   // Step (10)
 186   // d1 = a - S1 * S1 in f9
 187   (p6) fnma.s1 f9=f7,f7,f8
 188   nop.i 0;;
 189 } { .mfb
 190   nop.m 0
 191   // Step (11)
 192   // S = S1 + d1 * H1 in f7
 193   (p6) fma.s.s0 f8=f9,f12,f7
 194   (p6) br.ret.sptk   b0 ;;
 195 // END SINGLE PRECISION MINIMUM LATENCY SQUARE ROOT ALGORITHM
 196 } { .mfb
 197   nop.m 0
 198        mov   f8 = f7
 199   (p8) br.ret.sptk b0 ;;
 200 }
 201 //
 202 // This branch includes all those special values that are not negative,
 203 // with the result equal to frcpa(x)
 204 //
 205 GLOBAL_IEEE754_END(sqrtf)
 206
 207
 208 LOCAL_LIBM_ENTRY(__libm_error_region)
 209 .prologue
 210 { .mii
 211         add   GR_Parameter_Y=-32,sp             // Parameter 2 value
 212         mov   GR_Parameter_TAG = 50
 213 .save   ar.pfs,GR_SAVE_PFS
 214         mov  GR_SAVE_PFS=ar.pfs                 // Save ar.pfs
 215 }
 216 { .mfi
 217 .fframe 64
 218         add sp=-64,sp                           // Create new stack
 219         nop.f 0
 220         mov GR_SAVE_GP=gp                       // Save gp
 221 };;
 222 { .mmi
 223         stfs [GR_Parameter_Y] = FR_Y,16         // Store Parameter 2 on stack
 224         add GR_Parameter_X = 16,sp              // Parameter 1 address
 225 .save   b0, GR_SAVE_B0
 226         mov GR_SAVE_B0=b0                       // Save b0
 227 };;
 228 .body
 229 { .mib
 230         stfs [GR_Parameter_X] = FR_X            // Store Parameter 1 on stack
 231         add   GR_Parameter_RESULT = 0,GR_Parameter_Y
 232         nop.b 0                                 // Parameter 3 address
 233 }
 234 { .mib
 235         stfs [GR_Parameter_Y] = FR_RESULT       // Store Parameter 3 on stack
 236         add   GR_Parameter_Y = -16,GR_Parameter_Y
 237         br.call.sptk b0=__libm_error_support#   // Call error handling function
 238 };;
 239 { .mmi
 240         nop.m 0
 241         nop.m 0
 242         add   GR_Parameter_RESULT = 48,sp
 243 };;
 244 { .mmi
 245         ldfs  f8 = [GR_Parameter_RESULT]       // Get return result off stack
 246 .restore sp
 247         add   sp = 64,sp                       // Restore stack pointer
 248         mov   b0 = GR_SAVE_B0                  // Restore return address
 249 };;
 250 { .mib
 251         mov   gp = GR_SAVE_GP                  // Restore gp
 252         mov   ar.pfs = GR_SAVE_PFS             // Restore ar.pfs
 253         br.ret.sptk     b0                     // Return
 254 };;
 255
 256 LOCAL_LIBM_END(__libm_error_region)
 257
 258
 259 .type   __libm_error_support#,@function
 260 .global __libm_error_support#