softfloat/f32_sqrt.c

   1
   2 #include <stdbool.h>
   3 #include <stdint.h>
   4 #include "platform.h"
   5 #include "primitives.h"
   6 #include "internals.h"
   7 #include "specialize.h"
   8 #include "softfloat.h"
   9
  10 float32_t f32_sqrt( float32_t a )
  11 {
  12     union ui32_f32 uA;
  13     uint_fast32_t uiA;
  14     bool signA;
  15     int_fast16_t expA;
  16     uint_fast32_t sigA, uiZ;
  17     struct exp16_sig32 normExpSig;
  18     int_fast16_t expZ;
  19     uint_fast32_t sigZ;
  20     uint_fast64_t term, rem;
  21     union ui32_f32 uZ;
  22
  23     uA.f = a;
  24     uiA = uA.ui;
  25     signA = signF32UI( uiA );
  26     expA = expF32UI( uiA );
  27     sigA = fracF32UI( uiA );
  28     if ( expA == 0xFF ) {
  29         if ( sigA ) {
  30             uiZ = softfloat_propagateNaNF32UI( uiA, 0 );
  31             goto uiZ;
  32         }
  33         if ( ! signA ) return a;
  34         goto invalid;
  35     }
  36     if ( signA ) {
  37         if ( ! ( expA | sigA ) ) return a;
  38         goto invalid;
  39     }
  40     if ( ! expA ) {
  41         if ( ! sigA ) return a;
  42         normExpSig = softfloat_normSubnormalF32Sig( sigA );
  43         expA = normExpSig.exp;
  44         sigA = normExpSig.sig;
  45     }
  46     expZ = ( ( expA - 0x7F )>>1 ) + 0x7E;
  47     sigA = ( sigA | 0x00800000 )<<8;
  48     sigZ = softfloat_estimateSqrt32( expA, sigA ) + 2;
  49     if ( ( sigZ & 0x7F ) <= 5 ) {
  50         if ( sigZ < 2 ) {
  51             sigZ = 0x7FFFFFFF;
  52             goto roundPack;
  53         }
  54         sigA >>= expA & 1;
  55         term = (uint_fast64_t) sigZ * sigZ;
  56         rem = ( (uint_fast64_t) sigA<<32 ) - term;
  57         while ( UINT64_C( 0x8000000000000000 ) <= rem ) {
  58             --sigZ;
  59             rem += ( (uint_fast64_t) sigZ<<1 ) | 1;
  60         }
  61         sigZ |= ( rem != 0 );
  62     }
  63     sigZ = softfloat_shortShift32Right1Jam( sigZ );
  64  roundPack:
  65     return softfloat_roundPackToF32( 0, expZ, sigZ );
  66  invalid:
  67     softfloat_raiseFlags( softfloat_flag_invalid );
  68     uiZ = defaultNaNF32UI;
  69  uiZ:
  70     uZ.ui = uiZ;
  71     return uZ.f;
  72
  73 }
  74