From 95dbacb8fe7237b727789fb1fc1db16ff836ce87 Mon Sep 17 00:00:00 2001 From: Alessandro Andrioni Date: Sun, 19 May 2013 05:49:40 -0300 Subject: [PATCH 01/11] Add rounding modes to normal floating point arithmetic --- base/exports.jl | 8 ++++++++ base/rounding.jl | 44 ++++++++++++++++++++++++++++++++++++++++++++ base/sysimg.jl | 4 ++++ test/Makefile | 2 +- test/rounding.jl | 44 ++++++++++++++++++++++++++++++++++++++++++++ test/runtests.jl | 2 +- 6 files changed, 102 insertions(+), 2 deletions(-) create mode 100644 base/rounding.jl create mode 100644 test/rounding.jl diff --git a/base/exports.jl b/base/exports.jl index e05aec41870e6..12f20b86596e9 100644 --- a/base/exports.jl +++ b/base/exports.jl @@ -82,6 +82,11 @@ export Reverse, RevString, RopeString, + RoundDown, + RoundingMode, + RoundToNearest, + RoundToZero, + RoundUp, Schur, Set, SparseMatrixCSC, @@ -833,6 +838,9 @@ export get_bigfloat_rounding, set_bigfloat_rounding, with_bigfloat_rounding, + get_rounding, + set_rounding, + with_rounding, # statistics cor, diff --git a/base/rounding.jl b/base/rounding.jl new file mode 100644 index 0000000000000..be41ad90a0c35 --- /dev/null +++ b/base/rounding.jl @@ -0,0 +1,44 @@ +module Rounding + +export + RoundingMode, RoundToNearest, RoundToZero, RoundUp, RoundDown, + get_rounding, set_rounding, with_rounding + +## rounding modes ## +abstract RoundingMode +type RoundToNearest <: RoundingMode end +type RoundToZero <: RoundingMode end +type RoundUp <: RoundingMode end +type RoundDown <: RoundingMode end + +@unix_only set_rounding(::Type{RoundToNearest}) = ccall(:fesetround, Cint, (Cint, ), 0) +@unix_only set_rounding(::Type{RoundToZero}) = ccall(:fesetround, Cint, (Cint, ), 3072) +@unix_only set_rounding(::Type{RoundUp}) = ccall(:fesetround, Cint, (Cint, ), 2048) +@unix_only set_rounding(::Type{RoundDown}) = ccall(:fesetround, Cint, (Cint, ), 1024) + +@unix_only function get_rounding() + r = ccall(:fegetround, Cint, ()) + if r == 0 + return RoundToNearest + elseif r == 1024 + return RoundDown + elseif r == 2048 + return RoundUp + elseif r == 3072 + return RoundToZero + else + error() + end +end + +@unix_only function with_rounding{T<:RoundingMode}(f::Function, rounding::Type{T}) + old_rounding = get_rounding() + set_rounding(rounding) + try + return f() + finally + set_rounding(old_rounding) + end +end + +end #module \ No newline at end of file diff --git a/base/sysimg.jl b/base/sysimg.jl index 41d42d55000e0..af6b77bda2411 100644 --- a/base/sysimg.jl +++ b/base/sysimg.jl @@ -214,6 +214,10 @@ include("graphics.jl") include("profile.jl") importall .Profile +# rounding utilities +include("rounding.jl") +importall .Rounding + include = include_from_node1 # prime method cache with some things we know we'll need right after startup diff --git a/test/Makefile b/test/Makefile index 704190ee81e15..6b43042f9f98c 100644 --- a/test/Makefile +++ b/test/Makefile @@ -6,7 +6,7 @@ TESTS = all core keywordargs numbers strings unicode collections hashing \ math functional bigint sorting statistics spawn parallel arpack file \ git pkg resolve suitesparse complex version pollfd mpfr broadcast \ socket floatapprox priorityqueue readdlm regex float16 combinatorics \ - sysinfo + sysinfo rounding $(TESTS) :: @$(PRINT_JULIA) $(call spawn,$(JULIA_EXECUTABLE)) ./runtests.jl $@ diff --git a/test/rounding.jl b/test/rounding.jl new file mode 100644 index 0000000000000..e8956357651ab --- /dev/null +++ b/test/rounding.jl @@ -0,0 +1,44 @@ +# Small sanity tests to ensure changing the rounding of float functions work +using Base.Test +# a + b returns a number exactly between prevfloat(1.) and 1., so its +# final result depends strongly on the utilized rounding direction. +a = prevfloat(0.5) +b = 0.5 +c = 0x1p-54 +d = prevfloat(1.) + +# Default rounding direction, RoundToNearest +@test a + b === 1. +@test - a - b === -1. +@test a - b === -c +@test b - a === c + +# RoundToZero +with_rounding(RoundToZero) do + @test a + b === d + @test - a - b === -d + @test a - b === -c + @test b - a === c +end + +# Sanity check to see if we have returned to RoundToNearest +@test a + b === 1. +@test - a - b === -1. +@test a - b == -c +@test b - a == c + +# RoundUp +with_rounding(RoundUp) do + @test a + b === 1. + @test - a - b === -d + @test a - b === -c + @test b - a === c +end + +# RoundDown +with_rounding(RoundDown) do + @test a + b === d + @test - a - b === -1. + @test a - b === -c + @test b - a === c +end diff --git a/test/runtests.jl b/test/runtests.jl index a2b270081890b..4ef5791992c9a 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -6,7 +6,7 @@ testnames = ["core", "keywordargs", "numbers", "strings", "unicode", "arpack", "file", "suitesparse", "version", "resolve", "pollfd", "mpfr", "broadcast", "complex", "socket", "floatapprox", "readdlm", "regex", "float16", - "combinatorics", "sysinfo"] + "combinatorics", "sysinfo", "rounding"] tests = ARGS==["all"] ? testnames : ARGS From 823b0175ac2c98f27ba9bc61d7eb13e1ff16b4ef Mon Sep 17 00:00:00 2001 From: Alessandro Andrioni Date: Sun, 19 May 2013 06:04:04 -0300 Subject: [PATCH 02/11] Add Float32 rounding mode tests --- test/rounding.jl | 45 +++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 45 insertions(+) diff --git a/test/rounding.jl b/test/rounding.jl index e8956357651ab..bd0c13eb0ce07 100644 --- a/test/rounding.jl +++ b/test/rounding.jl @@ -1,5 +1,7 @@ # Small sanity tests to ensure changing the rounding of float functions work using Base.Test + +## Float64 checks # a + b returns a number exactly between prevfloat(1.) and 1., so its # final result depends strongly on the utilized rounding direction. a = prevfloat(0.5) @@ -42,3 +44,46 @@ with_rounding(RoundDown) do @test a - b === -c @test b - a === c end + +## Float32 checks + +a32 = prevfloat(0.5f0) +b32 = 0.5f0 +c32 = (1.f0 - prevfloat(1.f0))/2 +d32 = prevfloat(1.0f0) + +# Default rounding direction, RoundToNearest +@test a32 + b32 === 1.0f0 +@test - a32 - b32 === -1.0f0 +@test a32 - b32 === -c32 +@test b32 - a32 === c32 + +# RoundToZero +with_rounding(RoundToZero) do + @test a32 + b32 === d32 + @test - a32 - b32 === -d32 + @test a32 - b32 === -c32 + @test b32 - a32 === c32 +end + +# Sanity check to see if we have returned to RoundToNearest +@test a32 + b32 === 1.0f0 +@test - a32 - b32 === -1.0f0 +@test a32 - b32 == -c32 +@test b32 - a32 == c32 + +# RoundUp +with_rounding(RoundUp) do + @test a32 + b32 === 1.0f0 + @test - a32 - b32 === -d32 + @test a32 - b32 === -c32 + @test b32 - a32 === c32 +end + +# RoundDown +with_rounding(RoundDown) do + @test a32 + b32 === d32 + @test - a32 - b32 === -1.0f0 + @test a32 - b32 === -c32 + @test b32 - a32 === c32 +end From 1a132391a918a1dd8edb8de261f4ef26958d6a1d Mon Sep 17 00:00:00 2001 From: Alessandro Andrioni Date: Sun, 19 May 2013 06:48:06 -0300 Subject: [PATCH 03/11] Integrate the new rounding modes with BigFloats This changes the old rounding mode sketch to a more intuitive one, and using types instead of constants seems to be more logical. --- base/mpfr.jl | 41 ++++++++++++++++++++++++++--------------- base/sysimg.jl | 4 ++++ test/mpfr.jl | 6 +++--- 3 files changed, 33 insertions(+), 18 deletions(-) diff --git a/base/mpfr.jl b/base/mpfr.jl index 1ec78a5c7a7ea..e4373ad415659 100644 --- a/base/mpfr.jl +++ b/base/mpfr.jl @@ -19,7 +19,9 @@ import gamma, lgamma, digamma, erf, erfc, zeta, log1p, airyai, iceil, ifloor, itrunc, eps, signbit, sin, cos, tan, sec, csc, cot, acos, asin, atan, cosh, sinh, tanh, sech, csch, coth, acosh, asinh, atanh, atan2, - serialize, deserialize, inf, nan, hash + serialize, deserialize, inf, nan, hash, + RoundingMode, RoundDown, RoundingMode, RoundToNearest, RoundToZero, + RoundUp import Base.Math.lgamma_r @@ -27,11 +29,7 @@ const ROUNDING_MODE = [0] const DEFAULT_PRECISION = [256] # Rounding modes -const RoundToNearest = 0 -const RoundToZero = 1 -const RoundUp = 2 -const RoundDown = 3 -const RoundAwayZero = 4 +type RoundAwayZero <: RoundingMode end # Basic type and initialization definitions @@ -98,7 +96,7 @@ for to in (Int8, Int16, Int32, Int64) function convert(::Type{$to}, x::BigFloat) (isinteger(x) && (typemin($to) <= x <= typemax($to))) || throw(InexactError()) convert($to, ccall((:mpfr_get_si,:libmpfr), - Clong, (Ptr{BigFloat}, Int32), &x, RoundToZero)) + Clong, (Ptr{BigFloat}, Int32), &x, 0)) end end end @@ -108,7 +106,7 @@ for to in (Uint8, Uint16, Uint32, Uint64) function convert(::Type{$to}, x::BigFloat) (isinteger(x) && (typemin($to) <= x <= typemax($to))) || throw(InexactError()) convert($to, ccall((:mpfr_get_ui,:libmpfr), - Culong, (Ptr{BigFloat}, Int32), &x, RoundToZero)) + Culong, (Ptr{BigFloat}, Int32), &x, 0)) end end end @@ -597,13 +595,26 @@ function set_bigfloat_precision(x::Int) DEFAULT_PRECISION[end] = x end -get_bigfloat_rounding() = ROUNDING_MODE[end] -function set_bigfloat_rounding(x::Int) - if x < 0 || x > 4 - throw(DomainError()) +function get_bigfloat_rounding() + if ROUNDING_MODE[end] == 0 + return RoundToNearest + elseif ROUNDING_MODE[end] == 1 + return RoundToZero + elseif ROUNDING_MODE[end] == 2 + return RoundUp + elseif ROUNDING_MODE[end] == 3 + return RoundDown + elseif ROUNDING_MODE[end] == 4 + return RoundAwayZero + else + error("Invalid rounding mode") end - ROUNDING_MODE[end] = x end +set_bigfloat_rounding(::Type{RoundToNearest}) = ROUNDING_MODE[end] = 0 +set_bigfloat_rounding(::Type{RoundToZero}) = ROUNDING_MODE[end] = 1 +set_bigfloat_rounding(::Type{RoundUp}) = ROUNDING_MODE[end] = 2 +set_bigfloat_rounding(::Type{RoundDown}) = ROUNDING_MODE[end] = 3 +set_bigfloat_rounding(::Type{RoundAwayZero}) = ROUNDING_MODE[end] = 4 function copysign(x::BigFloat, y::BigFloat) z = BigFloat() @@ -635,7 +646,7 @@ end function itrunc(x::BigFloat) z = BigInt() - ccall((:mpfr_get_z, :libmpfr), Int32, (Ptr{BigInt}, Ptr{BigFloat}, Int32), &z, &x, RoundToZero) + ccall((:mpfr_get_z, :libmpfr), Int32, (Ptr{BigInt}, Ptr{BigFloat}, Int32), &z, &x, 0) return z end @@ -682,7 +693,7 @@ function with_bigfloat_precision(f::Function, precision::Integer) end end -function with_bigfloat_rounding(f::Function, rounding::Integer) +function with_bigfloat_rounding{T<:RoundingMode}(f::Function, rounding::Type{T}) old_rounding = get_bigfloat_rounding() set_bigfloat_rounding(rounding) try diff --git a/base/sysimg.jl b/base/sysimg.jl index af6b77bda2411..c05198d62ea8a 100644 --- a/base/sysimg.jl +++ b/base/sysimg.jl @@ -181,6 +181,10 @@ include("fftw.jl") include("dsp.jl") importall .DSP +# rounding utilities +include("rounding.jl") +importall .Rounding + # BigInts and BigFloats include("gmp.jl") importall .GMP diff --git a/test/mpfr.jl b/test/mpfr.jl index 257d09a7b85d8..cea88bb68e106 100644 --- a/test/mpfr.jl +++ b/test/mpfr.jl @@ -82,14 +82,14 @@ z = BigFloat(30) # rounding modes with_bigfloat_precision(4) do # default mode is round to nearest - down, up = with_bigfloat_rounding(MPFR.RoundToNearest) do + down, up = with_bigfloat_rounding(RoundToNearest) do BigFloat("0.0938"), BigFloat("0.102") end - with_bigfloat_rounding(MPFR.RoundDown) do + with_bigfloat_rounding(RoundDown) do @test BigFloat(0.1) == down @test BigFloat(0.1) != up end - with_bigfloat_rounding(MPFR.RoundUp) do + with_bigfloat_rounding(RoundUp) do @test BigFloat(0.1) != down @test BigFloat(0.1) == up end From 451f21e293f7d48c04283d4c36ada8f4c9094327 Mon Sep 17 00:00:00 2001 From: Alessandro Andrioni Date: Sun, 19 May 2013 06:57:07 -0300 Subject: [PATCH 04/11] Export RoundAwayZero too (even though it only works for BigFloats) --- base/exports.jl | 1 + base/mpfr.jl | 1 + 2 files changed, 2 insertions(+) diff --git a/base/exports.jl b/base/exports.jl index 12f20b86596e9..58d99b11ebbf2 100644 --- a/base/exports.jl +++ b/base/exports.jl @@ -82,6 +82,7 @@ export Reverse, RevString, RopeString, + RoundAwayZero, RoundDown, RoundingMode, RoundToNearest, diff --git a/base/mpfr.jl b/base/mpfr.jl index e4373ad415659..b5deabc37a865 100644 --- a/base/mpfr.jl +++ b/base/mpfr.jl @@ -2,6 +2,7 @@ module MPFR export BigFloat, + RoundAwayZero, get_bigfloat_precision, set_bigfloat_precision, with_bigfloat_precision, From 4f6ee0f1ac480ff14217f30f1549dbe8aa56a08e Mon Sep 17 00:00:00 2001 From: Alessandro Andrioni Date: Mon, 20 May 2013 05:50:05 -0300 Subject: [PATCH 05/11] Stop hardcoding the `fenv.h` constants As suggested by Viral in https://github.com/JuliaLang/julia/pull/3149 --- base/.gitignore | 1 + base/Makefile | 7 ++++++- base/rounding.jl | 21 +++++++++++---------- src/fenv_constants.h | 6 ++++++ 4 files changed, 24 insertions(+), 11 deletions(-) create mode 100644 src/fenv_constants.h diff --git a/base/.gitignore b/base/.gitignore index 77ecc7d193ca7..fd5e9b3ee9f4d 100644 --- a/base/.gitignore +++ b/base/.gitignore @@ -1,5 +1,6 @@ /pcre_h.jl /errno_h.jl /build_h.jl +/fenv_constants.jl /file_constants.jl /uv_constants.jl diff --git a/base/Makefile b/base/Makefile index f9a6bbd4d08ea..17114789c1ec2 100644 --- a/base/Makefile +++ b/base/Makefile @@ -3,6 +3,7 @@ include ../Make.inc PCRE_CONST = 0x[0-9a-fA-F]+|[-+]?\s*[0-9]+ +<<<<<<< HEAD # These are all the values needed for the RawVersionInfo struct version_string = $(shell cat ../VERSION) commit = $(shell git rev-parse HEAD 2>/dev/null) @@ -26,7 +27,7 @@ dirty = $(shell [ -z "$(shell git status --porcelain 2>/dev/null)" ] && echo "fa -all: pcre_h.jl errno_h.jl build_h.jl.phony file_constants.jl uv_constants.jl +all: pcre_h.jl errno_h.jl build_h.jl.phony fenv_constants.jl file_constants.jl uv_constants.jl pcre_h.jl: @$(PRINT_PERL) $(CPP) -dM $(shell $(PCRE_CONFIG) --prefix)/include/pcre.h | perl -nle '/^\s*#define\s+PCRE_(\w*)\s*\(?($(PCRE_CONST))\)?\s*$$/ and print "const $$1 = uint32($$2)"' | sort > $@ @@ -34,6 +35,9 @@ pcre_h.jl: errno_h.jl: @$(PRINT_PERL) echo '#include "errno.h"' | $(CPP) -dM - | perl -nle 'print "const $$1 = int32($$2)" if /^#define\s+(E\w+)\s+(\d+)\s*$$/' | sort > $@ +fenv_constants.jl: ../src/fenv_constants.h + $(QUIET_PERL) ${CC} -E -P -DJULIA ../src/fenv_constants.h | perl -nle 'print "$$1 0o$$2" if /^(\s*const\s+[A-z_]+\s+=)\s+(0[0-9]*)\s*$$/; print "$$1" if /^\s*(const\s+[A-z_]+\s+=\s+([1-9]|0x)[0-9A-z]*)\s*$$/' > $@ + file_constants.jl: ../src/file_constants.h @$(PRINT_PERL) $(CPP) -P -DJULIA ../src/file_constants.h | perl -nle 'print "$$1 0o$$2" if /^(\s*const\s+[A-z_]+\s+=)\s+(0[0-9]*)\s*$$/; print "$$1" if /^\s*(const\s+[A-z_]+\s+=\s+([1-9]|0x)[0-9A-z]*)\s*$$/' > $@ @@ -98,5 +102,6 @@ clean: rm -f pcre_h.jl rm -f errno_h.jl rm -f build_h.jl + rm -f fenv_constants.jl rm -f uv_constants.jl rm -f file_constants.jl diff --git a/base/rounding.jl b/base/rounding.jl index be41ad90a0c35..d75547ce05108 100644 --- a/base/rounding.jl +++ b/base/rounding.jl @@ -1,4 +1,5 @@ module Rounding +include("fenv_constants.jl") export RoundingMode, RoundToNearest, RoundToZero, RoundUp, RoundDown, @@ -11,27 +12,27 @@ type RoundToZero <: RoundingMode end type RoundUp <: RoundingMode end type RoundDown <: RoundingMode end -@unix_only set_rounding(::Type{RoundToNearest}) = ccall(:fesetround, Cint, (Cint, ), 0) -@unix_only set_rounding(::Type{RoundToZero}) = ccall(:fesetround, Cint, (Cint, ), 3072) -@unix_only set_rounding(::Type{RoundUp}) = ccall(:fesetround, Cint, (Cint, ), 2048) -@unix_only set_rounding(::Type{RoundDown}) = ccall(:fesetround, Cint, (Cint, ), 1024) +set_rounding(::Type{RoundToNearest}) = ccall(:fesetround, Cint, (Cint, ), JL_FE_TONEAREST) +set_rounding(::Type{RoundToZero}) = ccall(:fesetround, Cint, (Cint, ), JL_FE_TOWARDZERO) +set_rounding(::Type{RoundUp}) = ccall(:fesetround, Cint, (Cint, ), JL_FE_UPWARD) +set_rounding(::Type{RoundDown}) = ccall(:fesetround, Cint, (Cint, ), JL_FE_DOWNWARD) -@unix_only function get_rounding() +function get_rounding() r = ccall(:fegetround, Cint, ()) - if r == 0 + if r == JL_FE_TONEAREST return RoundToNearest - elseif r == 1024 + elseif r == JL_FE_DOWNWARD return RoundDown - elseif r == 2048 + elseif r == JL_FE_UPWARD return RoundUp - elseif r == 3072 + elseif r == JL_FE_TOWARDZERO return RoundToZero else error() end end -@unix_only function with_rounding{T<:RoundingMode}(f::Function, rounding::Type{T}) +function with_rounding{T<:RoundingMode}(f::Function, rounding::Type{T}) old_rounding = get_rounding() set_rounding(rounding) try diff --git a/src/fenv_constants.h b/src/fenv_constants.h new file mode 100644 index 0000000000000..64aaab5e18716 --- /dev/null +++ b/src/fenv_constants.h @@ -0,0 +1,6 @@ +#include + +const JL_FE_TONEAREST = FE_TONEAREST +const JL_FE_UPWARD = FE_UPWARD +const JL_FE_DOWNWARD = FE_DOWNWARD +const JL_FE_TOWARDZERO = FE_TOWARDZERO From 454f61a0523f98fa20f5d5f38ac1644e370c7f4b Mon Sep 17 00:00:00 2001 From: Alessandro Andrioni Date: Mon, 20 May 2013 06:21:22 -0300 Subject: [PATCH 06/11] Fix bug while parsing fenv_constants.h --- base/Makefile | 2 +- src/fenv_constants.h | 6 +++++- 2 files changed, 6 insertions(+), 2 deletions(-) diff --git a/base/Makefile b/base/Makefile index 17114789c1ec2..a53332cfd37e7 100644 --- a/base/Makefile +++ b/base/Makefile @@ -36,7 +36,7 @@ errno_h.jl: @$(PRINT_PERL) echo '#include "errno.h"' | $(CPP) -dM - | perl -nle 'print "const $$1 = int32($$2)" if /^#define\s+(E\w+)\s+(\d+)\s*$$/' | sort > $@ fenv_constants.jl: ../src/fenv_constants.h - $(QUIET_PERL) ${CC} -E -P -DJULIA ../src/fenv_constants.h | perl -nle 'print "$$1 0o$$2" if /^(\s*const\s+[A-z_]+\s+=)\s+(0[0-9]*)\s*$$/; print "$$1" if /^\s*(const\s+[A-z_]+\s+=\s+([1-9]|0x)[0-9A-z]*)\s*$$/' > $@ + $(QUIET_PERL) ${CC} -E -P -DJULIA ../src/fenv_constants.h | tail -n 8 > $@ file_constants.jl: ../src/file_constants.h @$(PRINT_PERL) $(CPP) -P -DJULIA ../src/file_constants.h | perl -nle 'print "$$1 0o$$2" if /^(\s*const\s+[A-z_]+\s+=)\s+(0[0-9]*)\s*$$/; print "$$1" if /^\s*(const\s+[A-z_]+\s+=\s+([1-9]|0x)[0-9A-z]*)\s*$$/' > $@ diff --git a/src/fenv_constants.h b/src/fenv_constants.h index 64aaab5e18716..9c9e96640b597 100644 --- a/src/fenv_constants.h +++ b/src/fenv_constants.h @@ -1,5 +1,9 @@ #include - +const JL_FE_INEXACT = FE_INEXACT +const JL_FE_UNDERFLOW = FE_UNDERFLOW +const JL_FE_OVERFLOW = FE_OVERFLOW +const JL_FE_DIVBYZERO = FE_DIVBYZERO +const JL_FE_INVALID = FE_INVALID const JL_FE_TONEAREST = FE_TONEAREST const JL_FE_UPWARD = FE_UPWARD const JL_FE_DOWNWARD = FE_DOWNWARD From 5a064b7778dd488ef76823fe21abffce082a0a3b Mon Sep 17 00:00:00 2001 From: Alessandro Andrioni Date: Mon, 26 Aug 2013 21:59:34 -0300 Subject: [PATCH 07/11] Renaming of the rounding modes --- base/exports.jl | 4 ++-- base/mpfr.jl | 14 +++++++------- base/rounding.jl | 8 ++++---- test/mpfr.jl | 2 +- test/rounding.jl | 8 ++++---- 5 files changed, 18 insertions(+), 18 deletions(-) diff --git a/base/exports.jl b/base/exports.jl index 58d99b11ebbf2..ae405cbabb087 100644 --- a/base/exports.jl +++ b/base/exports.jl @@ -82,10 +82,10 @@ export Reverse, RevString, RopeString, - RoundAwayZero, + RoundFromZero, RoundDown, RoundingMode, - RoundToNearest, + RoundNearest, RoundToZero, RoundUp, Schur, diff --git a/base/mpfr.jl b/base/mpfr.jl index b5deabc37a865..b7e2bca18a43c 100644 --- a/base/mpfr.jl +++ b/base/mpfr.jl @@ -2,7 +2,7 @@ module MPFR export BigFloat, - RoundAwayZero, + RoundFromZero, get_bigfloat_precision, set_bigfloat_precision, with_bigfloat_precision, @@ -21,7 +21,7 @@ import itrunc, eps, signbit, sin, cos, tan, sec, csc, cot, acos, asin, atan, cosh, sinh, tanh, sech, csch, coth, acosh, asinh, atanh, atan2, serialize, deserialize, inf, nan, hash, - RoundingMode, RoundDown, RoundingMode, RoundToNearest, RoundToZero, + RoundingMode, RoundDown, RoundingMode, RoundNearest, RoundToZero, RoundUp import Base.Math.lgamma_r @@ -30,7 +30,7 @@ const ROUNDING_MODE = [0] const DEFAULT_PRECISION = [256] # Rounding modes -type RoundAwayZero <: RoundingMode end +type RoundFromZero <: RoundingMode end # Basic type and initialization definitions @@ -598,7 +598,7 @@ end function get_bigfloat_rounding() if ROUNDING_MODE[end] == 0 - return RoundToNearest + return RoundNearest elseif ROUNDING_MODE[end] == 1 return RoundToZero elseif ROUNDING_MODE[end] == 2 @@ -606,16 +606,16 @@ function get_bigfloat_rounding() elseif ROUNDING_MODE[end] == 3 return RoundDown elseif ROUNDING_MODE[end] == 4 - return RoundAwayZero + return RoundFromZero else error("Invalid rounding mode") end end -set_bigfloat_rounding(::Type{RoundToNearest}) = ROUNDING_MODE[end] = 0 +set_bigfloat_rounding(::Type{RoundNearest}) = ROUNDING_MODE[end] = 0 set_bigfloat_rounding(::Type{RoundToZero}) = ROUNDING_MODE[end] = 1 set_bigfloat_rounding(::Type{RoundUp}) = ROUNDING_MODE[end] = 2 set_bigfloat_rounding(::Type{RoundDown}) = ROUNDING_MODE[end] = 3 -set_bigfloat_rounding(::Type{RoundAwayZero}) = ROUNDING_MODE[end] = 4 +set_bigfloat_rounding(::Type{RoundFromZero}) = ROUNDING_MODE[end] = 4 function copysign(x::BigFloat, y::BigFloat) z = BigFloat() diff --git a/base/rounding.jl b/base/rounding.jl index d75547ce05108..736a8bbeed06e 100644 --- a/base/rounding.jl +++ b/base/rounding.jl @@ -2,17 +2,17 @@ module Rounding include("fenv_constants.jl") export - RoundingMode, RoundToNearest, RoundToZero, RoundUp, RoundDown, + RoundingMode, RoundNearest, RoundToZero, RoundUp, RoundDown, get_rounding, set_rounding, with_rounding ## rounding modes ## abstract RoundingMode -type RoundToNearest <: RoundingMode end +type RoundNearest <: RoundingMode end type RoundToZero <: RoundingMode end type RoundUp <: RoundingMode end type RoundDown <: RoundingMode end -set_rounding(::Type{RoundToNearest}) = ccall(:fesetround, Cint, (Cint, ), JL_FE_TONEAREST) +set_rounding(::Type{RoundNearest}) = ccall(:fesetround, Cint, (Cint, ), JL_FE_TONEAREST) set_rounding(::Type{RoundToZero}) = ccall(:fesetround, Cint, (Cint, ), JL_FE_TOWARDZERO) set_rounding(::Type{RoundUp}) = ccall(:fesetround, Cint, (Cint, ), JL_FE_UPWARD) set_rounding(::Type{RoundDown}) = ccall(:fesetround, Cint, (Cint, ), JL_FE_DOWNWARD) @@ -20,7 +20,7 @@ set_rounding(::Type{RoundDown}) = ccall(:fesetround, Cint, (Cint, ), JL_FE_DOWNW function get_rounding() r = ccall(:fegetround, Cint, ()) if r == JL_FE_TONEAREST - return RoundToNearest + return RoundNearest elseif r == JL_FE_DOWNWARD return RoundDown elseif r == JL_FE_UPWARD diff --git a/test/mpfr.jl b/test/mpfr.jl index cea88bb68e106..18b71cf66df83 100644 --- a/test/mpfr.jl +++ b/test/mpfr.jl @@ -82,7 +82,7 @@ z = BigFloat(30) # rounding modes with_bigfloat_precision(4) do # default mode is round to nearest - down, up = with_bigfloat_rounding(RoundToNearest) do + down, up = with_bigfloat_rounding(RoundNearest) do BigFloat("0.0938"), BigFloat("0.102") end with_bigfloat_rounding(RoundDown) do diff --git a/test/rounding.jl b/test/rounding.jl index bd0c13eb0ce07..0194dd3567b36 100644 --- a/test/rounding.jl +++ b/test/rounding.jl @@ -9,7 +9,7 @@ b = 0.5 c = 0x1p-54 d = prevfloat(1.) -# Default rounding direction, RoundToNearest +# Default rounding direction, RoundNearest @test a + b === 1. @test - a - b === -1. @test a - b === -c @@ -23,7 +23,7 @@ with_rounding(RoundToZero) do @test b - a === c end -# Sanity check to see if we have returned to RoundToNearest +# Sanity check to see if we have returned to RoundNearest @test a + b === 1. @test - a - b === -1. @test a - b == -c @@ -52,7 +52,7 @@ b32 = 0.5f0 c32 = (1.f0 - prevfloat(1.f0))/2 d32 = prevfloat(1.0f0) -# Default rounding direction, RoundToNearest +# Default rounding direction, RoundNearest @test a32 + b32 === 1.0f0 @test - a32 - b32 === -1.0f0 @test a32 - b32 === -c32 @@ -66,7 +66,7 @@ with_rounding(RoundToZero) do @test b32 - a32 === c32 end -# Sanity check to see if we have returned to RoundToNearest +# Sanity check to see if we have returned to RoundNearest @test a32 + b32 === 1.0f0 @test - a32 - b32 === -1.0f0 @test a32 - b32 == -c32 From 54a4055b6dd7ed93531552d25c083a0a25033ba6 Mon Sep 17 00:00:00 2001 From: Alessandro Andrioni Date: Mon, 26 Aug 2013 22:09:53 -0300 Subject: [PATCH 08/11] Fixes post-rebase --- base/Makefile | 1 - base/fenv_constants.jl | 8 ++++++++ base/sysimg.jl | 4 ---- 3 files changed, 8 insertions(+), 5 deletions(-) create mode 100644 base/fenv_constants.jl diff --git a/base/Makefile b/base/Makefile index a53332cfd37e7..94cf333528b52 100644 --- a/base/Makefile +++ b/base/Makefile @@ -3,7 +3,6 @@ include ../Make.inc PCRE_CONST = 0x[0-9a-fA-F]+|[-+]?\s*[0-9]+ -<<<<<<< HEAD # These are all the values needed for the RawVersionInfo struct version_string = $(shell cat ../VERSION) commit = $(shell git rev-parse HEAD 2>/dev/null) diff --git a/base/fenv_constants.jl b/base/fenv_constants.jl new file mode 100644 index 0000000000000..82088b2ed321d --- /dev/null +++ b/base/fenv_constants.jl @@ -0,0 +1,8 @@ +const JL_FE_UNDERFLOW = 0x0010 +const JL_FE_OVERFLOW = 0x0008 +const JL_FE_DIVBYZERO = 0x0004 +const JL_FE_INVALID = 0x0001 +const JL_FE_TONEAREST = 0x0000 +const JL_FE_UPWARD = 0x0800 +const JL_FE_DOWNWARD = 0x0400 +const JL_FE_TOWARDZERO = 0x0c00 diff --git a/base/sysimg.jl b/base/sysimg.jl index c05198d62ea8a..e64691e0d1b88 100644 --- a/base/sysimg.jl +++ b/base/sysimg.jl @@ -218,10 +218,6 @@ include("graphics.jl") include("profile.jl") importall .Profile -# rounding utilities -include("rounding.jl") -importall .Rounding - include = include_from_node1 # prime method cache with some things we know we'll need right after startup From e314eba1f5e068b085ebc4b53e366c3b81a1f4fc Mon Sep 17 00:00:00 2001 From: Alessandro Andrioni Date: Tue, 27 Aug 2013 21:15:52 -0300 Subject: [PATCH 09/11] Remove fenv_constants.jl, as it is auto-generated --- base/fenv_constants.jl | 8 -------- 1 file changed, 8 deletions(-) delete mode 100644 base/fenv_constants.jl diff --git a/base/fenv_constants.jl b/base/fenv_constants.jl deleted file mode 100644 index 82088b2ed321d..0000000000000 --- a/base/fenv_constants.jl +++ /dev/null @@ -1,8 +0,0 @@ -const JL_FE_UNDERFLOW = 0x0010 -const JL_FE_OVERFLOW = 0x0008 -const JL_FE_DIVBYZERO = 0x0004 -const JL_FE_INVALID = 0x0001 -const JL_FE_TONEAREST = 0x0000 -const JL_FE_UPWARD = 0x0800 -const JL_FE_DOWNWARD = 0x0400 -const JL_FE_TOWARDZERO = 0x0c00 From 285b2389d80de9b0b778dfd7db192bb43fd66415 Mon Sep 17 00:00:00 2001 From: Alessandro Andrioni Date: Tue, 27 Aug 2013 22:06:22 -0300 Subject: [PATCH 10/11] Documentation for the rounding modes (manual + stdlib) It also includes some updates to some related sections about BigFloats. --- .../integers-and-floating-point-numbers.rst | 45 ++++++++++++++++++- doc/stdlib/base.rst | 21 ++++++++- 2 files changed, 63 insertions(+), 3 deletions(-) diff --git a/doc/manual/integers-and-floating-point-numbers.rst b/doc/manual/integers-and-floating-point-numbers.rst index 5a666da88ed35..dda4f56c4a2fc 100644 --- a/doc/manual/integers-and-floating-point-numbers.rst +++ b/doc/manual/integers-and-floating-point-numbers.rst @@ -447,6 +447,26 @@ argument respectively: :: This example highlights the general principle that the adjacent representable floating-point numbers also have adjacent binary integer representations. +Rounding modes +~~~~~~~~~~~~~~ + +If a number doesn't have an exact floating-point representation, it must be +rounded to an appropriate representable value, however, if wanted, the manner +in which this rounding is done can be changed according to the rounding modes +presented in the `IEEE 754 standard `_:: + + + julia> 1.1 + 0.1 + 1.2000000000000002 + + julia> with_rounding(RoundDown) do + 1.1 + 0.1 + end + 1.2 + +The default mode used is always ``RoundNearest``, which rounds to the nearest +representable value, with ties rounded towards the nearest value with an even +least significant bit. Background and References ~~~~~~~~~~~~~~~~~~~~~~~~~ @@ -503,10 +523,10 @@ Once created, they participate in arithmetic with all other numeric types thanks 123456789012345678901234567891 julia> BigFloat("1.23456789012345678901") - 1.23456789012345678901 + 1.234567890123456789010000000000000000000000000000000000000000000000000000000004e+00 with 256 bits of precision julia> BigFloat(2.0^66) / 3 - 24595658764946068821.3 + 2.459565876494606882133333333333333333333333333333333333333333333333333333333344e+19 with 256 bits of precision julia> factorial(BigInt(40)) 815915283247897734345611269596115894272000000000 @@ -531,6 +551,27 @@ However, type promotion between the primitive types above and julia> typeof(y) BigInt + +The default precision (in number of bits of the significand) and rounding +mode of `BigFloat` operations can be changed, and all further calculations +will take these changes in account:: + + julia> with_bigfloat_rounding(RoundUp) do + BigFloat(1) + BigFloat("0.1") + end + 1.100000000000000000000000000000000000000000000000000000000000000000000000000003e+00 with 256 bits of precision + + julia> with_bigfloat_rounding(RoundDown) do + BigFloat(1) + BigFloat("0.1") + end + 1.099999999999999999999999999999999999999999999999999999999999999999999999999986e+00 with 256 bits of precision + + julia> with_bigfloat_precision(40) do + BigFloat(1) + BigFloat("0.1") + end + 1.0999999999985e+00 with 40 bits of precision + + .. _man-numeric-literal-coefficients: diff --git a/doc/stdlib/base.rst b/doc/stdlib/base.rst index c6c5f31a6504a..a7bf5fa65a7c1 100644 --- a/doc/stdlib/base.rst +++ b/doc/stdlib/base.rst @@ -2902,6 +2902,25 @@ Numbers Create an arbitrary precision floating point number. ``x`` may be an ``Integer``, a ``Float64``, a ``String`` or a ``BigInt``. The usual mathematical operators are defined for this type, and results are promoted to a ``BigFloat``. +.. function:: get_rounding() + + Get the current floating point rounding mode. Valid modes are ``RoundNearest``, ``RoundToZero``, ``RoundUp`` and ``RoundDown``. + +.. function:: set_rounding(mode) + + Set the floating point rounding mode. See ``get_rounding`` for available modes + +.. function:: with_rounding(f::Function,mode) + + Change the floating point rounding mode for the duration of ``f``. It is logically equivalent to:: + + old = get_rounding() + set_rounding(mode) + f() + set_rounding(old) + + See ``get_rounding`` for available rounding modes. + Integers ~~~~~~~~ @@ -2990,7 +3009,7 @@ The `BigFloat` type implements arbitrary-precision floating-point aritmetic usin .. function:: get_bigfloat_rounding() - Get the current BigFloat rounding mode. Valid modes are ``RoundToNearest``, ``RoundToZero``, ``RoundUp``, ``RoundDown``, ``RoundAwayZero`` + Get the current BigFloat rounding mode. Valid modes are ``RoundNearest``, ``RoundToZero``, ``RoundUp``, ``RoundDown``, ``RoundFromZero`` .. function:: set_bigfloat_rounding(mode) From 773b60e442472828c76c33209c3e7b0581588a44 Mon Sep 17 00:00:00 2001 From: Alessandro Andrioni Date: Wed, 28 Aug 2013 11:39:51 -0300 Subject: [PATCH 11/11] Fix issue with fenv.h and Ubuntu 12.04 If the values from fenv.h aren't found by the C preprocessor, it uses the default ones for x86/x64 --- base/Makefile | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/base/Makefile b/base/Makefile index 94cf333528b52..328fff691dd5f 100644 --- a/base/Makefile +++ b/base/Makefile @@ -35,7 +35,7 @@ errno_h.jl: @$(PRINT_PERL) echo '#include "errno.h"' | $(CPP) -dM - | perl -nle 'print "const $$1 = int32($$2)" if /^#define\s+(E\w+)\s+(\d+)\s*$$/' | sort > $@ fenv_constants.jl: ../src/fenv_constants.h - $(QUIET_PERL) ${CC} -E -P -DJULIA ../src/fenv_constants.h | tail -n 8 > $@ + @$(PRINT_PERL) ${CPP} -P -DJULIA ../src/fenv_constants.h | tail -n 8 | perl -ple 's/\sFE_UN\w+/ 0x10/g; s/\sFE_O\w+/ 0x08/g; s/\sFE_DI\w+/ 0x04/g; s/\sFE_INV\w+/ 0x01/g; s/\sFE_TON\w+/ 0x00/g; s/\sFE_UP\w+/ 0x800/g; s/\sFE_DO\w+/ 0x400/g; s/\sFE_TOW\w+/ 0xc00/g' > $@ file_constants.jl: ../src/file_constants.h @$(PRINT_PERL) $(CPP) -P -DJULIA ../src/file_constants.h | perl -nle 'print "$$1 0o$$2" if /^(\s*const\s+[A-z_]+\s+=)\s+(0[0-9]*)\s*$$/; print "$$1" if /^\s*(const\s+[A-z_]+\s+=\s+([1-9]|0x)[0-9A-z]*)\s*$$/' > $@