libgo/go/math/remainder.go

   1 // Copyright 2010 The Go Authors. All rights reserved.
   2 // Use of this source code is governed by a BSD-style
   3 // license that can be found in the LICENSE file.
   4
   5 package math
   6
   7 // The original C code and the comment below are from
   8 // FreeBSD's /usr/src/lib/msun/src/e_remainder.c and came
   9 // with this notice. The go code is a simplified version of
  10 // the original C.
  11 //
  12 // ====================================================
  13 // Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
  14 //
  15 // Developed at SunPro, a Sun Microsystems, Inc. business.
  16 // Permission to use, copy, modify, and distribute this
  17 // software is freely granted, provided that this notice
  18 // is preserved.
  19 // ====================================================
  20 //
  21 // __ieee754_remainder(x,y)
  22 // Return :
  23 //      returns  x REM y  =  x - [x/y]*y  as if in infinite
  24 //      precision arithmetic, where [x/y] is the (infinite bit)
  25 //      integer nearest x/y (in half way cases, choose the even one).
  26 // Method :
  27 //      Based on Mod() returning  x - [x/y]chopped * y  exactly.
  28
  29 // Remainder returns the IEEE 754 floating-point remainder of x/y.
  30 //
  31 // Special cases are:
  32 //      Remainder(±Inf, y) = NaN
  33 //      Remainder(NaN, y) = NaN
  34 //      Remainder(x, 0) = NaN
  35 //      Remainder(x, ±Inf) = x
  36 //      Remainder(x, NaN) = NaN
  37 func Remainder(x, y float64) float64 {
  38         if haveArchRemainder {
  39                 return archRemainder(x, y)
  40         }
  41         return remainder(x, y)
  42 }
  43
  44 func remainder(x, y float64) float64 {
  45         const (
  46                 Tiny    = 4.45014771701440276618e-308 // 0x0020000000000000
  47                 HalfMax = MaxFloat64 / 2
  48         )
  49         // special cases
  50         switch {
  51         case IsNaN(x) || IsNaN(y) || IsInf(x, 0) || y == 0:
  52                 return NaN()
  53         case IsInf(y, 0):
  54                 return x
  55         }
  56         sign := false
  57         if x < 0 {
  58                 x = -x
  59                 sign = true
  60         }
  61         if y < 0 {
  62                 y = -y
  63         }
  64         if x == y {
  65                 if sign {
  66                         zero := 0.0
  67                         return -zero
  68                 }
  69                 return 0
  70         }
  71         if y <= HalfMax {
  72                 x = Mod(x, y+y) // now x < 2y
  73         }
  74         if y < Tiny {
  75                 if x+x > y {
  76                         x -= y
  77                         if x+x >= y {
  78                                 x -= y
  79                         }
  80                 }
  81         } else {
  82                 yHalf := 0.5 * y
  83                 if x > yHalf {
  84                         x -= y
  85                         if x >= yHalf {
  86                                 x -= y
  87                         }
  88                 }
  89         }
  90         if sign {
  91                 x = -x
  92         }
  93         return x
  94 }