From 4ee0544a157090ddd087b36109d292cd174bae7c Mon Sep 17 00:00:00 2001 From: xleroy Date: Fri, 2 Aug 2013 08:05:18 +0000 Subject: Merge of Flocq version 2.2.0. More precise modeling of NaNs. git-svn-id: https://yquem.inria.fr/compcert/svn/compcert/trunk@2303 fca1b0fc-160b-0410-b1d3-a4f43f01ea2e --- .depend | 38 +- Changelog | 2 + Makefile | 22 +- arm/Nan.v | 34 ++ common/Values.v | 5 - extraction/extraction.v | 10 +- flocq/Appli/Fappli_IEEE.v | 412 ++++++++++------ flocq/Appli/Fappli_IEEE_bits.v | 68 ++- flocq/Appli/Fappli_rnd_odd.v | 979 ++++++++++++++++++++++++++++++++++++++ flocq/Calc/Fcalc_bracket.v | 4 +- flocq/Calc/Fcalc_digits.v | 4 +- flocq/Calc/Fcalc_div.v | 4 +- flocq/Calc/Fcalc_ops.v | 4 +- flocq/Calc/Fcalc_round.v | 4 +- flocq/Calc/Fcalc_sqrt.v | 4 +- flocq/Core/Fcore.v | 4 +- flocq/Core/Fcore_FIX.v | 4 +- flocq/Core/Fcore_FLT.v | 4 +- flocq/Core/Fcore_FLX.v | 4 +- flocq/Core/Fcore_FTZ.v | 4 +- flocq/Core/Fcore_Raux.v | 16 +- flocq/Core/Fcore_Zaux.v | 4 +- flocq/Core/Fcore_defs.v | 4 +- flocq/Core/Fcore_digits.v | 4 +- flocq/Core/Fcore_float_prop.v | 4 +- flocq/Core/Fcore_generic_fmt.v | 103 +++- flocq/Core/Fcore_rnd.v | 4 +- flocq/Core/Fcore_rnd_ne.v | 23 +- flocq/Core/Fcore_ulp.v | 278 ++++++++++- flocq/Flocq_version.v | 6 +- flocq/Prop/Fprop_Sterbenz.v | 4 +- flocq/Prop/Fprop_div_sqrt_error.v | 4 +- flocq/Prop/Fprop_mult_error.v | 4 +- flocq/Prop/Fprop_plus_error.v | 4 +- flocq/Prop/Fprop_relative.v | 166 ++++++- ia32/Nan.v | 28 ++ lib/Floats.v | 444 +++++++++++++---- powerpc/Nan.v | 28 ++ 38 files changed, 2378 insertions(+), 364 deletions(-) create mode 100644 arm/Nan.v create mode 100644 flocq/Appli/Fappli_rnd_odd.v create mode 100644 ia32/Nan.v create mode 100644 powerpc/Nan.v diff --git a/.depend b/.depend index ca707a4..0e6c556 100644 --- a/.depend +++ b/.depend @@ -8,6 +8,7 @@ lib/Ordered.vo lib/Ordered.glob lib/Ordered.v.beautified: lib/Ordered.v lib/Coql lib/Iteration.vo lib/Iteration.glob lib/Iteration.v.beautified: lib/Iteration.v lib/Axioms.vo lib/Coqlib.vo lib/Wfsimpl.vo lib/Integers.vo lib/Integers.glob lib/Integers.v.beautified: lib/Integers.v lib/Coqlib.vo lib/Floats.vo lib/Floats.glob lib/Floats.v.beautified: lib/Floats.v lib/Axioms.vo lib/Coqlib.vo lib/Integers.vo flocq/Appli/Fappli_IEEE.vo flocq/Appli/Fappli_IEEE_bits.vo flocq/Core/Fcore.vo flocq/Calc/Fcalc_round.vo flocq/Calc/Fcalc_bracket.vo flocq/Prop/Fprop_Sterbenz.vo +$(ARCH)/Nan.vo $(ARCH)/Nan.glob $(ARCH)/Nan.v.beautified: $(ARCH)/Nan.v flocq/Appli/Fappli_IEEE.vo flocq/Appli/Fappli_IEEE_bits.vo lib/Floats.vo lib/Integers.vo lib/Parmov.vo lib/Parmov.glob lib/Parmov.v.beautified: lib/Parmov.v lib/Axioms.vo lib/Coqlib.vo lib/UnionFind.vo lib/UnionFind.glob lib/UnionFind.v.beautified: lib/UnionFind.v lib/Coqlib.vo lib/Wfsimpl.vo lib/Wfsimpl.glob lib/Wfsimpl.v.beautified: lib/Wfsimpl.v lib/Axioms.vo @@ -107,31 +108,32 @@ cfrontend/Cminorgen.vo cfrontend/Cminorgen.glob cfrontend/Cminorgen.v.beautified cfrontend/Cminorgenproof.vo cfrontend/Cminorgenproof.glob cfrontend/Cminorgenproof.v.beautified: cfrontend/Cminorgenproof.v lib/Coqlib.vo lib/Intv.vo common/Errors.vo lib/Maps.vo common/AST.vo lib/Integers.vo lib/Floats.vo common/Values.vo common/Memory.vo common/Events.vo common/Globalenvs.vo common/Smallstep.vo common/Switch.vo cfrontend/Csharpminor.vo backend/Cminor.vo cfrontend/Cminorgen.vo driver/Compiler.vo driver/Compiler.glob driver/Compiler.v.beautified: driver/Compiler.v lib/Coqlib.vo common/Errors.vo common/AST.vo common/Smallstep.vo cfrontend/Csyntax.vo cfrontend/Csem.vo cfrontend/Cstrategy.vo cfrontend/Cexec.vo cfrontend/Clight.vo cfrontend/Csharpminor.vo backend/Cminor.vo backend/CminorSel.vo backend/RTL.vo backend/LTL.vo backend/Linear.vo backend/Mach.vo $(ARCH)/Asm.vo cfrontend/Initializers.vo cfrontend/SimplExpr.vo cfrontend/SimplLocals.vo cfrontend/Cshmgen.vo cfrontend/Cminorgen.vo backend/Selection.vo backend/RTLgen.vo backend/Tailcall.vo backend/Inlining.vo backend/Renumber.vo backend/Constprop.vo backend/CSE.vo backend/Allocation.vo backend/Tunneling.vo backend/Linearize.vo backend/CleanupLabels.vo backend/Stacking.vo $(ARCH)/Asmgen.vo cfrontend/SimplExprproof.vo cfrontend/SimplLocalsproof.vo cfrontend/Cshmgenproof.vo cfrontend/Cminorgenproof.vo backend/Selectionproof.vo backend/RTLgenproof.vo backend/Tailcallproof.vo backend/Inliningproof.vo backend/Renumberproof.vo backend/Constpropproof.vo backend/CSEproof.vo backend/Allocproof.vo backend/Tunnelingproof.vo backend/Linearizeproof.vo backend/CleanupLabelsproof.vo backend/Stackingproof.vo $(ARCH)/Asmgenproof.vo driver/Complements.vo driver/Complements.glob driver/Complements.v.beautified: driver/Complements.v lib/Coqlib.vo common/AST.vo lib/Integers.vo common/Values.vo common/Events.vo common/Globalenvs.vo common/Smallstep.vo common/Behaviors.vo cfrontend/Csyntax.vo cfrontend/Csem.vo cfrontend/Cstrategy.vo cfrontend/Clight.vo backend/Cminor.vo backend/RTL.vo $(ARCH)/Asm.vo driver/Compiler.vo common/Errors.vo -flocq/Core/Fcore_float_prop.vo flocq/Core/Fcore_float_prop.glob flocq/Core/Fcore_float_prop.v.beautified: flocq/Core/Fcore_float_prop.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo +flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_Raux.glob flocq/Core/Fcore_Raux.v.beautified: flocq/Core/Fcore_Raux.v flocq/Core/Fcore_Zaux.vo flocq/Core/Fcore_Zaux.vo flocq/Core/Fcore_Zaux.glob flocq/Core/Fcore_Zaux.v.beautified: flocq/Core/Fcore_Zaux.v -flocq/Core/Fcore_rnd_ne.vo flocq/Core/Fcore_rnd_ne.glob flocq/Core/Fcore_rnd_ne.v.beautified: flocq/Core/Fcore_rnd_ne.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_rnd.vo flocq/Core/Fcore_generic_fmt.vo flocq/Core/Fcore_float_prop.vo flocq/Core/Fcore_ulp.vo -flocq/Core/Fcore_FTZ.vo flocq/Core/Fcore_FTZ.glob flocq/Core/Fcore_FTZ.v.beautified: flocq/Core/Fcore_FTZ.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_rnd.vo flocq/Core/Fcore_generic_fmt.vo flocq/Core/Fcore_float_prop.vo flocq/Core/Fcore_FLX.vo -flocq/Core/Fcore_FLT.vo flocq/Core/Fcore_FLT.glob flocq/Core/Fcore_FLT.v.beautified: flocq/Core/Fcore_FLT.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_rnd.vo flocq/Core/Fcore_generic_fmt.vo flocq/Core/Fcore_float_prop.vo flocq/Core/Fcore_FLX.vo flocq/Core/Fcore_FIX.vo flocq/Core/Fcore_rnd_ne.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_defs.glob flocq/Core/Fcore_defs.v.beautified: flocq/Core/Fcore_defs.v flocq/Core/Fcore_Raux.vo -flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_Raux.glob flocq/Core/Fcore_Raux.v.beautified: flocq/Core/Fcore_Raux.v flocq/Core/Fcore_Zaux.vo -flocq/Core/Fcore_ulp.vo flocq/Core/Fcore_ulp.glob flocq/Core/Fcore_ulp.v.beautified: flocq/Core/Fcore_ulp.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_rnd.vo flocq/Core/Fcore_generic_fmt.vo flocq/Core/Fcore_float_prop.vo -flocq/Core/Fcore_rnd.vo flocq/Core/Fcore_rnd.glob flocq/Core/Fcore_rnd.v.beautified: flocq/Core/Fcore_rnd.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo -flocq/Core/Fcore_FLX.vo flocq/Core/Fcore_FLX.glob flocq/Core/Fcore_FLX.v.beautified: flocq/Core/Fcore_FLX.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_rnd.vo flocq/Core/Fcore_generic_fmt.vo flocq/Core/Fcore_float_prop.vo flocq/Core/Fcore_FIX.vo flocq/Core/Fcore_rnd_ne.vo -flocq/Core/Fcore_FIX.vo flocq/Core/Fcore_FIX.glob flocq/Core/Fcore_FIX.v.beautified: flocq/Core/Fcore_FIX.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_rnd.vo flocq/Core/Fcore_generic_fmt.vo flocq/Core/Fcore_rnd_ne.vo flocq/Core/Fcore_digits.vo flocq/Core/Fcore_digits.glob flocq/Core/Fcore_digits.v.beautified: flocq/Core/Fcore_digits.v flocq/Core/Fcore_Zaux.vo +flocq/Core/Fcore_float_prop.vo flocq/Core/Fcore_float_prop.glob flocq/Core/Fcore_float_prop.v.beautified: flocq/Core/Fcore_float_prop.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo +flocq/Core/Fcore_FIX.vo flocq/Core/Fcore_FIX.glob flocq/Core/Fcore_FIX.v.beautified: flocq/Core/Fcore_FIX.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_rnd.vo flocq/Core/Fcore_generic_fmt.vo flocq/Core/Fcore_rnd_ne.vo +flocq/Core/Fcore_FLT.vo flocq/Core/Fcore_FLT.glob flocq/Core/Fcore_FLT.v.beautified: flocq/Core/Fcore_FLT.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_rnd.vo flocq/Core/Fcore_generic_fmt.vo flocq/Core/Fcore_float_prop.vo flocq/Core/Fcore_FLX.vo flocq/Core/Fcore_FIX.vo flocq/Core/Fcore_rnd_ne.vo +flocq/Core/Fcore_FLX.vo flocq/Core/Fcore_FLX.glob flocq/Core/Fcore_FLX.v.beautified: flocq/Core/Fcore_FLX.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_rnd.vo flocq/Core/Fcore_generic_fmt.vo flocq/Core/Fcore_float_prop.vo flocq/Core/Fcore_FIX.vo flocq/Core/Fcore_rnd_ne.vo +flocq/Core/Fcore_FTZ.vo flocq/Core/Fcore_FTZ.glob flocq/Core/Fcore_FTZ.v.beautified: flocq/Core/Fcore_FTZ.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_rnd.vo flocq/Core/Fcore_generic_fmt.vo flocq/Core/Fcore_float_prop.vo flocq/Core/Fcore_FLX.vo flocq/Core/Fcore_generic_fmt.vo flocq/Core/Fcore_generic_fmt.glob flocq/Core/Fcore_generic_fmt.v.beautified: flocq/Core/Fcore_generic_fmt.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_rnd.vo flocq/Core/Fcore_float_prop.vo +flocq/Core/Fcore_rnd.vo flocq/Core/Fcore_rnd.glob flocq/Core/Fcore_rnd.v.beautified: flocq/Core/Fcore_rnd.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo +flocq/Core/Fcore_rnd_ne.vo flocq/Core/Fcore_rnd_ne.glob flocq/Core/Fcore_rnd_ne.v.beautified: flocq/Core/Fcore_rnd_ne.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_rnd.vo flocq/Core/Fcore_generic_fmt.vo flocq/Core/Fcore_float_prop.vo flocq/Core/Fcore_ulp.vo +flocq/Core/Fcore_ulp.vo flocq/Core/Fcore_ulp.glob flocq/Core/Fcore_ulp.v.beautified: flocq/Core/Fcore_ulp.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_rnd.vo flocq/Core/Fcore_generic_fmt.vo flocq/Core/Fcore_float_prop.vo flocq/Core/Fcore.vo flocq/Core/Fcore.glob flocq/Core/Fcore.v.beautified: flocq/Core/Fcore.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_float_prop.vo flocq/Core/Fcore_rnd.vo flocq/Core/Fcore_generic_fmt.vo flocq/Core/Fcore_rnd_ne.vo flocq/Core/Fcore_FIX.vo flocq/Core/Fcore_FLX.vo flocq/Core/Fcore_FLT.vo flocq/Core/Fcore_ulp.vo -flocq/Prop/Fprop_Sterbenz.vo flocq/Prop/Fprop_Sterbenz.glob flocq/Prop/Fprop_Sterbenz.v.beautified: flocq/Prop/Fprop_Sterbenz.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_generic_fmt.vo flocq/Calc/Fcalc_ops.vo -flocq/Prop/Fprop_mult_error.vo flocq/Prop/Fprop_mult_error.glob flocq/Prop/Fprop_mult_error.v.beautified: flocq/Prop/Fprop_mult_error.v flocq/Core/Fcore.vo flocq/Calc/Fcalc_ops.vo -flocq/Prop/Fprop_relative.vo flocq/Prop/Fprop_relative.glob flocq/Prop/Fprop_relative.v.beautified: flocq/Prop/Fprop_relative.v flocq/Core/Fcore.vo -flocq/Prop/Fprop_div_sqrt_error.vo flocq/Prop/Fprop_div_sqrt_error.glob flocq/Prop/Fprop_div_sqrt_error.v.beautified: flocq/Prop/Fprop_div_sqrt_error.v flocq/Core/Fcore.vo flocq/Calc/Fcalc_ops.vo flocq/Prop/Fprop_relative.vo -flocq/Prop/Fprop_plus_error.vo flocq/Prop/Fprop_plus_error.glob flocq/Prop/Fprop_plus_error.v.beautified: flocq/Prop/Fprop_plus_error.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_float_prop.vo flocq/Core/Fcore_generic_fmt.vo flocq/Calc/Fcalc_ops.vo -flocq/Calc/Fcalc_ops.vo flocq/Calc/Fcalc_ops.glob flocq/Calc/Fcalc_ops.v.beautified: flocq/Calc/Fcalc_ops.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_float_prop.vo flocq/Calc/Fcalc_bracket.vo flocq/Calc/Fcalc_bracket.glob flocq/Calc/Fcalc_bracket.v.beautified: flocq/Calc/Fcalc_bracket.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_float_prop.vo -flocq/Calc/Fcalc_sqrt.vo flocq/Calc/Fcalc_sqrt.glob flocq/Calc/Fcalc_sqrt.v.beautified: flocq/Calc/Fcalc_sqrt.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_digits.vo flocq/Core/Fcore_float_prop.vo flocq/Calc/Fcalc_bracket.vo flocq/Calc/Fcalc_digits.vo +flocq/Calc/Fcalc_digits.vo flocq/Calc/Fcalc_digits.glob flocq/Calc/Fcalc_digits.v.beautified: flocq/Calc/Fcalc_digits.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_float_prop.vo flocq/Core/Fcore_digits.vo flocq/Calc/Fcalc_div.vo flocq/Calc/Fcalc_div.glob flocq/Calc/Fcalc_div.v.beautified: flocq/Calc/Fcalc_div.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_float_prop.vo flocq/Core/Fcore_digits.vo flocq/Calc/Fcalc_bracket.vo flocq/Calc/Fcalc_digits.vo +flocq/Calc/Fcalc_ops.vo flocq/Calc/Fcalc_ops.glob flocq/Calc/Fcalc_ops.v.beautified: flocq/Calc/Fcalc_ops.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_float_prop.vo flocq/Calc/Fcalc_round.vo flocq/Calc/Fcalc_round.glob flocq/Calc/Fcalc_round.v.beautified: flocq/Calc/Fcalc_round.v flocq/Core/Fcore.vo flocq/Core/Fcore_digits.vo flocq/Calc/Fcalc_bracket.vo flocq/Calc/Fcalc_digits.vo -flocq/Calc/Fcalc_digits.vo flocq/Calc/Fcalc_digits.glob flocq/Calc/Fcalc_digits.v.beautified: flocq/Calc/Fcalc_digits.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_float_prop.vo flocq/Core/Fcore_digits.vo -flocq/Appli/Fappli_IEEE_bits.vo flocq/Appli/Fappli_IEEE_bits.glob flocq/Appli/Fappli_IEEE_bits.v.beautified: flocq/Appli/Fappli_IEEE_bits.v flocq/Core/Fcore.vo flocq/Core/Fcore_digits.vo flocq/Calc/Fcalc_digits.vo flocq/Appli/Fappli_IEEE.vo +flocq/Calc/Fcalc_sqrt.vo flocq/Calc/Fcalc_sqrt.glob flocq/Calc/Fcalc_sqrt.v.beautified: flocq/Calc/Fcalc_sqrt.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_digits.vo flocq/Core/Fcore_float_prop.vo flocq/Calc/Fcalc_bracket.vo flocq/Calc/Fcalc_digits.vo +flocq/Prop/Fprop_div_sqrt_error.vo flocq/Prop/Fprop_div_sqrt_error.glob flocq/Prop/Fprop_div_sqrt_error.v.beautified: flocq/Prop/Fprop_div_sqrt_error.v flocq/Core/Fcore.vo flocq/Calc/Fcalc_ops.vo flocq/Prop/Fprop_relative.vo +flocq/Prop/Fprop_mult_error.vo flocq/Prop/Fprop_mult_error.glob flocq/Prop/Fprop_mult_error.v.beautified: flocq/Prop/Fprop_mult_error.v flocq/Core/Fcore.vo flocq/Calc/Fcalc_ops.vo +flocq/Prop/Fprop_plus_error.vo flocq/Prop/Fprop_plus_error.glob flocq/Prop/Fprop_plus_error.v.beautified: flocq/Prop/Fprop_plus_error.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_float_prop.vo flocq/Core/Fcore_generic_fmt.vo flocq/Calc/Fcalc_ops.vo +flocq/Prop/Fprop_relative.vo flocq/Prop/Fprop_relative.glob flocq/Prop/Fprop_relative.v.beautified: flocq/Prop/Fprop_relative.v flocq/Core/Fcore.vo +flocq/Prop/Fprop_Sterbenz.vo flocq/Prop/Fprop_Sterbenz.glob flocq/Prop/Fprop_Sterbenz.v.beautified: flocq/Prop/Fprop_Sterbenz.v flocq/Core/Fcore_Raux.vo flocq/Core/Fcore_defs.vo flocq/Core/Fcore_generic_fmt.vo flocq/Calc/Fcalc_ops.vo +flocq/Appli/Fappli_rnd_odd.vo flocq/Appli/Fappli_rnd_odd.glob flocq/Appli/Fappli_rnd_odd.v.beautified: flocq/Appli/Fappli_rnd_odd.v flocq/Core/Fcore.vo flocq/Calc/Fcalc_ops.vo flocq/Appli/Fappli_IEEE.vo flocq/Appli/Fappli_IEEE.glob flocq/Appli/Fappli_IEEE.v.beautified: flocq/Appli/Fappli_IEEE.v flocq/Core/Fcore.vo flocq/Core/Fcore_digits.vo flocq/Calc/Fcalc_digits.vo flocq/Calc/Fcalc_round.vo flocq/Calc/Fcalc_bracket.vo flocq/Calc/Fcalc_ops.vo flocq/Calc/Fcalc_div.vo flocq/Calc/Fcalc_sqrt.vo flocq/Prop/Fprop_relative.vo +flocq/Appli/Fappli_IEEE_bits.vo flocq/Appli/Fappli_IEEE_bits.glob flocq/Appli/Fappli_IEEE_bits.v.beautified: flocq/Appli/Fappli_IEEE_bits.v flocq/Core/Fcore.vo flocq/Core/Fcore_digits.vo flocq/Calc/Fcalc_digits.vo flocq/Appli/Fappli_IEEE.vo exportclight/Clightdefs.vo exportclight/Clightdefs.glob exportclight/Clightdefs.v.beautified: exportclight/Clightdefs.v lib/Integers.vo lib/Floats.vo common/AST.vo cfrontend/Ctypes.vo cfrontend/Cop.vo cfrontend/Clight.vo diff --git a/Changelog b/Changelog index 657cb03..b7e37a7 100644 --- a/Changelog +++ b/Changelog @@ -1,6 +1,8 @@ Development trunk: ================== +- More precise modeling of not-a-numbers (NaNs) in floating-point + arithmetic. - Optimize integer divisions by positive constants, turning them into multiply-high and shifts. - Avoid double rounding issues in conversion from 64-bit integers diff --git a/Makefile b/Makefile index 1694457..66765d6 100644 --- a/Makefile +++ b/Makefile @@ -50,20 +50,22 @@ VPATH=$(DIRS) GPATH=$(DIRS) # Flocq -FLOCQ_CORE=Fcore_float_prop.v Fcore_Zaux.v Fcore_rnd_ne.v Fcore_FTZ.v \ - Fcore_FLT.v Fcore_defs.v Fcore_Raux.v Fcore_ulp.v Fcore_rnd.v Fcore_FLX.v \ - Fcore_FIX.v Fcore_digits.v Fcore_generic_fmt.v Fcore.v -FLOCQ_PROP=Fprop_Sterbenz.v Fprop_mult_error.v Fprop_relative.v \ - Fprop_div_sqrt_error.v Fprop_plus_error.v -FLOCQ_CALC=Fcalc_ops.v Fcalc_bracket.v Fcalc_sqrt.v Fcalc_div.v Fcalc_round.v \ - Fcalc_digits.v -FLOCQ_APPLI=Fappli_IEEE_bits.v Fappli_IEEE.v -FLOCQ=$(FLOCQ_CORE) $(FLOCQ_PROP) $(FLOCQ_CALC) $(FLOCQ_APPLI) + +FLOCQ=\ + Fcore_Raux.v Fcore_Zaux.v Fcore_defs.v Fcore_digits.v \ + Fcore_float_prop.v Fcore_FIX.v Fcore_FLT.v Fcore_FLX.v \ + Fcore_FTZ.v Fcore_generic_fmt.v Fcore_rnd.v Fcore_rnd_ne.v \ + Fcore_ulp.v Fcore.v \ + Fcalc_bracket.v Fcalc_digits.v Fcalc_div.v Fcalc_ops.v \ + Fcalc_round.v Fcalc_sqrt.v \ + Fprop_div_sqrt_error.v Fprop_mult_error.v Fprop_plus_error.v \ + Fprop_relative.v Fprop_Sterbenz.v \ + Fappli_rnd_odd.v Fappli_IEEE.v Fappli_IEEE_bits.v # General-purpose libraries (in lib/) LIB=Axioms.v Coqlib.v Intv.v Maps.v Heaps.v Lattice.v Ordered.v \ - Iteration.v Integers.v Floats.v Parmov.v UnionFind.v Wfsimpl.v \ + Iteration.v Integers.v Floats.v Nan.v Parmov.v UnionFind.v Wfsimpl.v \ Postorder.v FSetAVLplus.v # Parts common to the front-ends and the back-end (in common/) diff --git a/arm/Nan.v b/arm/Nan.v new file mode 100644 index 0000000..d66d6ec --- /dev/null +++ b/arm/Nan.v @@ -0,0 +1,34 @@ +Require Import Fappli_IEEE. +Require Import Fappli_IEEE_bits. +Require Import Floats. +Require Import ZArith. +Require Import Integers. + +(* Needed to break a circular reference after extraction *) +Definition transform_quiet_pl := + Eval unfold Float.transform_quiet_pl in Float.transform_quiet_pl. + +Program Definition default_pl : bool * nan_pl 53 := (false, nat_iter 51 xO xH). + +Definition binop_pl (pl1 pl2:binary64) : bool*nan_pl 53 := + match pl1, pl2 with + | B754_nan s1 pl1, B754_nan s2 pl2 => + if (Pos.testbit (proj1_sig pl1) 51 && (* pl2 is sNan but pl1 is qNan *) + negb (Pos.testbit (proj1_sig pl2) 51))%bool + then (s2, transform_quiet_pl pl2) + else (s1, transform_quiet_pl pl1) + | B754_nan s1 pl1, _ => (s1, transform_quiet_pl pl1) + | _, B754_nan s2 pl2 => (s2, transform_quiet_pl pl2) + | _, _ => default_pl + end. + +Theorem binop_propagate1: Float.binop_propagate1_prop binop_pl. +Proof. + repeat intro. destruct f1, f2; try discriminate; simpl; + destruct Pos.testbit; reflexivity. +Qed. + +Theorem binop_propagate2: Float.binop_propagate2_prop binop_pl. +Proof. + repeat intro. destruct f1, f2, f3; try discriminate; simpl; reflexivity. +Qed. diff --git a/common/Values.v b/common/Values.v index 670f785..b9594fc 100644 --- a/common/Values.v +++ b/common/Values.v @@ -1051,11 +1051,6 @@ Proof. intros; destruct x; simpl; auto. decEq. apply Int.rolm_zero. Qed. -Theorem addf_commut: forall x y, addf x y = addf y x. -Proof. - destruct x; destruct y; simpl; auto. decEq. apply Float.addf_commut. -Qed. - Theorem negate_cmp_bool: forall c x y, cmp_bool (negate_comparison c) x y = option_map negb (cmp_bool c x y). Proof. diff --git a/extraction/extraction.v b/extraction/extraction.v index 6bc0234..ba1ca58 100644 --- a/extraction/extraction.v +++ b/extraction/extraction.v @@ -11,6 +11,7 @@ (* *********************************************************************) Require Wfsimpl. +Require Nan. Require AST. Require Iteration. Require Floats. @@ -30,6 +31,10 @@ Require Import ExtrOcamlString. (* Wfsimpl *) Extraction Inline Wfsimpl.Fix Wfsimpl.Fixm. +(* Floats *) +Extract Constant Floats.Float.binop_pl => + "Nan.binop_pl". + (* AST *) Extract Constant AST.ident_of_string => "fun s -> Camlcoq.intern_string (Camlcoq.camlstring_of_coqstring s)". @@ -129,6 +134,5 @@ Separate Extraction Conventions1.dummy_int_reg Conventions1.dummy_float_reg RTL.instr_defs RTL.instr_uses Machregs.mregs_for_operation Machregs.mregs_for_builtin - Machregs.two_address_op. - - + Machregs.two_address_op + Nan.binop_pl. diff --git a/flocq/Appli/Fappli_IEEE.v b/flocq/Appli/Fappli_IEEE.v index 63b150f..5e9897f 100644 --- a/flocq/Appli/Fappli_IEEE.v +++ b/flocq/Appli/Fappli_IEEE.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2010-2011 Sylvie Boldo +Copyright (C) 2010-2013 Sylvie Boldo #
# -Copyright (C) 2010-2011 Guillaume Melquiond +Copyright (C) 2010-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public @@ -33,7 +33,7 @@ Section AnyRadix. Inductive full_float := | F754_zero : bool -> full_float | F754_infinity : bool -> full_float - | F754_nan : full_float + | F754_nan : bool -> positive -> full_float | F754_finite : bool -> positive -> Z -> full_float. Definition FF2R beta x := @@ -67,16 +67,20 @@ Definition bounded m e := Definition valid_binary x := match x with | F754_finite _ m e => bounded m e + | F754_nan _ pl => (Z_of_nat' (S (digits2_Pnat pl)) true end. (** Basic type used for representing binary FP numbers. Note that there is exactly one such object per FP datum. NaNs do not have any payload. They cannot be distinguished. *) + +Definition nan_pl := {pl | (Z_of_nat' (S (digits2_Pnat pl)) binary_float | B754_infinity : bool -> binary_float - | B754_nan : binary_float + | B754_nan : bool -> nan_pl -> binary_float | B754_finite : bool -> forall (m : positive) (e : Z), bounded m e = true -> binary_float. @@ -85,7 +89,7 @@ Definition FF2B x := | F754_finite s m e => B754_finite s m e | F754_infinity s => fun _ => B754_infinity s | F754_zero s => fun _ => B754_zero s - | F754_nan => fun _ => B754_nan + | F754_nan b pl => fun H => B754_nan b (exist _ pl H) end. Definition B2FF x := @@ -93,7 +97,7 @@ Definition B2FF x := | B754_finite s m e _ => F754_finite s m e | B754_infinity s => F754_infinity s | B754_zero s => F754_zero s - | B754_nan => F754_nan + | B754_nan b (exist pl _) => F754_nan b pl end. Definition radix2 := Build_radix 2 (refl_equal true). @@ -108,30 +112,30 @@ Theorem FF2R_B2FF : forall x, FF2R radix2 (B2FF x) = B2R x. Proof. -now intros [sx|sx| |sx mx ex Hx]. +now intros [sx|sx|sx [plx Hplx]|sx mx ex Hx]. Qed. Theorem B2FF_FF2B : forall x Hx, B2FF (FF2B x Hx) = x. Proof. -now intros [sx|sx| |sx mx ex] Hx. +now intros [sx|sx|sx plx|sx mx ex] Hx. Qed. Theorem valid_binary_B2FF : forall x, valid_binary (B2FF x) = true. Proof. -now intros [sx|sx| |sx mx ex Hx]. +now intros [sx|sx|sx [plx Hplx]|sx mx ex Hx]. Qed. Theorem FF2B_B2FF : forall x H, FF2B (B2FF x) H = x. Proof. -intros [sx|sx| |sx mx ex Hx] H ; try easy. -apply f_equal. -apply eqbool_irrelevance. +intros [sx|sx|sx [plx Hplx]|sx mx ex Hx] H ; try easy. +simpl. apply f_equal, f_equal, eqbool_irrelevance. +apply f_equal, eqbool_irrelevance. Qed. Theorem FF2B_B2FF_valid : @@ -146,7 +150,7 @@ Theorem B2R_FF2B : forall x Hx, B2R (FF2B x Hx) = FF2R radix2 x. Proof. -now intros [sx|sx| |sx mx ex] Hx. +now intros [sx|sx|sx plx|sx mx ex] Hx. Qed. Theorem match_FF2B : @@ -154,17 +158,17 @@ Theorem match_FF2B : match FF2B x Hx return T with | B754_zero sx => fz sx | B754_infinity sx => fi sx - | B754_nan => fn + | B754_nan b (exist p _) => fn b p | B754_finite sx mx ex _ => ff sx mx ex end = match x with | F754_zero sx => fz sx | F754_infinity sx => fi sx - | F754_nan => fn + | F754_nan b p => fn b p | F754_finite sx mx ex => ff sx mx ex end. Proof. -now intros T fz fi fn ff [sx|sx| |sx mx ex] Hx. +now intros T fz fi fn ff [sx|sx|sx plx|sx mx ex] Hx. Qed. Theorem canonic_canonic_mantissa : @@ -189,19 +193,28 @@ Theorem generic_format_B2R : forall x, generic_format radix2 fexp (B2R x). Proof. -intros [sx|sx| |sx mx ex Hx] ; try apply generic_format_0. +intros [sx|sx|sx plx|sx mx ex Hx] ; try apply generic_format_0. simpl. apply generic_format_canonic. apply canonic_canonic_mantissa. now destruct (andb_prop _ _ Hx) as (H, _). Qed. +Theorem FLT_format_B2R : + forall x, + FLT_format radix2 emin prec (B2R x). +Proof with auto with typeclass_instances. +intros x. +apply FLT_format_generic... +apply generic_format_B2R. +Qed. + Theorem B2FF_inj : forall x y : binary_float, B2FF x = B2FF y -> x = y. Proof. -intros [sx|sx| |sx mx ex Hx] [sy|sy| |sy my ey Hy] ; try easy. +intros [sx|sx|sx [plx Hplx]|sx mx ex Hx] [sy|sy|sy [ply Hply]|sy my ey Hy] ; try easy. (* *) intros H. now inversion H. @@ -212,11 +225,18 @@ now inversion H. intros H. inversion H. clear H. +revert Hplx. +rewrite H2. +intros Hx. +apply f_equal, f_equal, eqbool_irrelevance. +(* *) +intros H. +inversion H. +clear H. revert Hx. rewrite H2, H3. intros Hx. -apply f_equal. -apply eqbool_irrelevance. +apply f_equal, eqbool_irrelevance. Qed. Definition is_finite_strict f := @@ -265,6 +285,29 @@ apply f_equal. apply eqbool_irrelevance. Qed. +Definition Bsign x := + match x with + | B754_nan s _ => s + | B754_zero s => s + | B754_infinity s => s + | B754_finite s _ _ _ => s + end. + +Definition sign_FF x := + match x with + | F754_nan s _ => s + | F754_zero s => s + | F754_infinity s => s + | F754_finite s _ _ => s + end. + +Theorem Bsign_FF2B : + forall x H, + Bsign (FF2B x H) = sign_FF x. +Proof. +now intros [sx|sx|sx plx|sx mx ex] H. +Qed. + Definition is_finite f := match f with | B754_finite _ _ _ _ => true @@ -290,42 +333,93 @@ Theorem is_finite_FF_B2FF : forall x, is_finite_FF (B2FF x) = is_finite x. Proof. +now intros [| |? []|]. +Qed. + +Theorem B2R_Bsign_inj: + forall x y : binary_float, + is_finite x = true -> + is_finite y = true -> + B2R x = B2R y -> + Bsign x = Bsign y -> + x = y. +Proof. +intros. destruct x, y; try (apply B2R_inj; now eauto). +- simpl in H2. congruence. +- symmetry in H1. apply Rmult_integral in H1. + destruct H1. apply eq_Z2R with (n:=0%Z) in H1. destruct b0; discriminate H1. + simpl in H1. pose proof (bpow_gt_0 radix2 e). + rewrite H1 in H3. apply Rlt_irrefl in H3. destruct H3. +- apply Rmult_integral in H1. + destruct H1. apply eq_Z2R with (n:=0%Z) in H1. destruct b; discriminate H1. + simpl in H1. pose proof (bpow_gt_0 radix2 e). + rewrite H1 in H3. apply Rlt_irrefl in H3. destruct H3. +Qed. + +Definition is_nan f := + match f with + | B754_nan _ _ => true + | _ => false + end. + +Definition is_nan_FF f := + match f with + | F754_nan _ _ => true + | _ => false + end. + +Theorem is_nan_FF2B : + forall x Hx, + is_nan (FF2B x Hx) = is_nan_FF x. +Proof. now intros [| | |]. Qed. -Definition Bopp x := +Theorem is_nan_FF_B2FF : + forall x, + is_nan_FF (B2FF x) = is_nan x. +Proof. +now intros [| |? []|]. +Qed. + +Definition Bopp opp_nan x := match x with - | B754_nan => x + | B754_nan sx plx => + let '(sres, plres) := opp_nan sx plx in B754_nan sres plres | B754_infinity sx => B754_infinity (negb sx) | B754_finite sx mx ex Hx => B754_finite (negb sx) mx ex Hx | B754_zero sx => B754_zero (negb sx) end. Theorem Bopp_involutive : - forall x, Bopp (Bopp x) = x. + forall opp_nan x, + is_nan x = false -> + Bopp opp_nan (Bopp opp_nan x) = x. Proof. -now intros [sx|sx| |sx mx ex Hx] ; simpl ; try rewrite Bool.negb_involutive. +now intros opp_nan [sx|sx|sx plx|sx mx ex Hx] ; simpl ; try rewrite Bool.negb_involutive. Qed. Theorem B2R_Bopp : - forall x, - B2R (Bopp x) = (- B2R x)%R. + forall opp_nan x, + B2R (Bopp opp_nan x) = (- B2R x)%R. Proof. -intros [sx|sx| |sx mx ex Hx] ; apply sym_eq ; try apply Ropp_0. +intros opp_nan [sx|sx|sx plx|sx mx ex Hx]; apply sym_eq ; try apply Ropp_0. +simpl. destruct opp_nan. apply Ropp_0. simpl. rewrite <- F2R_opp. now case sx. Qed. - -Theorem is_finite_Bopp: forall x, - is_finite (Bopp x) = is_finite x. +Theorem is_finite_Bopp : + forall opp_nan x, + is_finite (Bopp opp_nan x) = is_finite x. Proof. -now intros [| | |]. +intros opp_nan [| | |] ; try easy. +intros s pl. +simpl. +now case opp_nan. Qed. - - Theorem bounded_lt_emax : forall mx ex, bounded mx ex = true -> @@ -361,7 +455,7 @@ Theorem abs_B2R_lt_emax : forall x, (Rabs (B2R x) < bpow radix2 emax)%R. Proof. -intros [sx|sx| |sx mx ex Hx] ; simpl ; try ( rewrite Rabs_R0 ; apply bpow_gt_0 ). +intros [sx|sx|sx plx|sx mx ex Hx] ; simpl ; try ( rewrite Rabs_R0 ; apply bpow_gt_0 ). rewrite <- F2R_Zabs, abs_cond_Zopp. now apply bounded_lt_emax. Qed. @@ -638,7 +732,7 @@ Definition binary_round_aux mode sx mx ex lx := match shr_m mrs'' with | Z0 => F754_zero sx | Zpos m => if Zle_bool e'' (emax - prec) then F754_finite sx m e'' else binary_overflow mode sx - | _ => F754_nan (* dummy *) + | _ => F754_nan false xH (* dummy *) end. Theorem binary_round_aux_correct : @@ -649,7 +743,7 @@ Theorem binary_round_aux_correct : valid_binary z = true /\ if Rlt_bool (Rabs (round radix2 fexp (round_mode mode) x)) (bpow radix2 emax) then FF2R radix2 z = round radix2 fexp (round_mode mode) x /\ - is_finite_FF z = true + is_finite_FF z = true /\ sign_FF z = Rlt_bool x 0 else z = binary_overflow mode (Rlt_bool x 0). Proof with auto with typeclass_instances. @@ -818,29 +912,6 @@ exact inbetween_int_UP_sign. exact inbetween_int_NA_sign. Qed. -Definition Bsign x := - match x with - | B754_nan => false - | B754_zero s => s - | B754_infinity s => s - | B754_finite s _ _ _ => s - end. - -Definition sign_FF x := - match x with - | F754_nan => false - | F754_zero s => s - | F754_infinity s => s - | F754_finite s _ _ => s - end. - -Theorem Bsign_FF2B : - forall x H, - Bsign (FF2B x H) = sign_FF x. -Proof. -now intros [sx|sx| |sx mx ex] H. -Qed. - (** Multiplication *) Lemma Bmult_correct_aux : @@ -851,7 +922,7 @@ Lemma Bmult_correct_aux : valid_binary z = true /\ if Rlt_bool (Rabs (round radix2 fexp (round_mode m) (x * y))) (bpow radix2 emax) then FF2R radix2 z = round radix2 fexp (round_mode m) (x * y) /\ - is_finite_FF z = true + is_finite_FF z = true /\ sign_FF z = xorb sx sy else z = binary_overflow m (xorb sx sy). Proof. @@ -896,15 +967,15 @@ apply Rlt_bool_false. now apply F2R_ge_0_compat. Qed. -Definition Bmult m x y := +Definition Bmult mult_nan m x y := + let f pl := B754_nan (fst pl) (snd pl) in match x, y with - | B754_nan, _ => x - | _, B754_nan => y + | B754_nan _ _, _ | _, B754_nan _ _ => f (mult_nan x y) | B754_infinity sx, B754_infinity sy => B754_infinity (xorb sx sy) | B754_infinity sx, B754_finite sy _ _ _ => B754_infinity (xorb sx sy) | B754_finite sx _ _ _, B754_infinity sy => B754_infinity (xorb sx sy) - | B754_infinity _, B754_zero _ => B754_nan - | B754_zero _, B754_infinity _ => B754_nan + | B754_infinity _, B754_zero _ => f (mult_nan x y) + | B754_zero _, B754_infinity _ => f (mult_nan x y) | B754_finite sx _ _ _, B754_zero sy => B754_zero (xorb sx sy) | B754_zero sx, B754_finite sy _ _ _ => B754_zero (xorb sx sy) | B754_zero sx, B754_zero sy => B754_zero (xorb sx sy) @@ -913,36 +984,40 @@ Definition Bmult m x y := end. Theorem Bmult_correct : - forall m x y, + forall mult_nan m x y, if Rlt_bool (Rabs (round radix2 fexp (round_mode m) (B2R x * B2R y))) (bpow radix2 emax) then - B2R (Bmult m x y) = round radix2 fexp (round_mode m) (B2R x * B2R y) /\ - is_finite (Bmult m x y) = andb (is_finite x) (is_finite y) + B2R (Bmult mult_nan m x y) = round radix2 fexp (round_mode m) (B2R x * B2R y) /\ + is_finite (Bmult mult_nan m x y) = andb (is_finite x) (is_finite y) /\ + (is_nan (Bmult mult_nan m x y) = false -> + Bsign (Bmult mult_nan m x y) = xorb (Bsign x) (Bsign y)) else - B2FF (Bmult m x y) = binary_overflow m (xorb (Bsign x) (Bsign y)). + B2FF (Bmult mult_nan m x y) = binary_overflow m (xorb (Bsign x) (Bsign y)). Proof. -intros m [sx|sx| |sx mx ex Hx] [sy|sy| |sy my ey Hy] ; - try ( rewrite ?Rmult_0_r, ?Rmult_0_l, round_0, Rabs_R0, Rlt_bool_true ; [ split ; apply refl_equal | apply bpow_gt_0 | auto with typeclass_instances ] ). +intros mult_nan m [sx|sx|sx plx|sx mx ex Hx] [sy|sy|sy ply|sy my ey Hy] ; + try ( rewrite ?Rmult_0_r, ?Rmult_0_l, round_0, Rabs_R0, Rlt_bool_true ; [ now repeat constructor | apply bpow_gt_0 | now auto with typeclass_instances ] ). simpl. case Bmult_correct_aux. intros H1. case Rlt_bool. -intros (H2, H3). +intros (H2, (H3, H4)). split. now rewrite B2R_FF2B. +split. now rewrite is_finite_FF2B. +rewrite Bsign_FF2B. auto. intros H2. now rewrite B2FF_FF2B. Qed. -Definition Bmult_FF m x y := +Definition Bmult_FF mult_nan m x y := + let f pl := F754_nan (fst pl) (snd pl) in match x, y with - | F754_nan, _ => x - | _, F754_nan => y + | F754_nan _ _, _ | _, F754_nan _ _ => f (mult_nan x y) | F754_infinity sx, F754_infinity sy => F754_infinity (xorb sx sy) | F754_infinity sx, F754_finite sy _ _ => F754_infinity (xorb sx sy) | F754_finite sx _ _, F754_infinity sy => F754_infinity (xorb sx sy) - | F754_infinity _, F754_zero _ => F754_nan - | F754_zero _, F754_infinity _ => F754_nan + | F754_infinity _, F754_zero _ => f (mult_nan x y) + | F754_zero _, F754_infinity _ => f (mult_nan x y) | F754_finite sx _ _, F754_zero sy => F754_zero (xorb sx sy) | F754_zero sx, F754_finite sy _ _ => F754_zero (xorb sx sy) | F754_zero sx, F754_zero sy => F754_zero (xorb sx sy) @@ -951,14 +1026,20 @@ Definition Bmult_FF m x y := end. Theorem B2FF_Bmult : + forall mult_nan mult_nan_ff, forall m x y, - B2FF (Bmult m x y) = Bmult_FF m (B2FF x) (B2FF y). + mult_nan_ff (B2FF x) (B2FF y) = (let '(sr, exist plr _) := mult_nan x y in (sr, plr)) -> + B2FF (Bmult mult_nan m x y) = Bmult_FF mult_nan_ff m (B2FF x) (B2FF y). Proof. -intros m [sx|sx| |sx mx ex Hx] [sy|sy| |sy my ey Hy] ; try easy. +intros mult_nan mult_nan_ff m x y Hmult_nan. +unfold Bmult_FF. rewrite Hmult_nan. +destruct x as [sx|sx|sx [plx Hplx]|sx mx ex Hx], y as [sy|sy|sy [ply Hply]|sy my ey Hy] ; + simpl; try match goal with |- context [mult_nan ?x ?y] => + destruct (mult_nan x y) as [? []] end; + try easy. apply B2FF_FF2B. Qed. - Definition shl_align mx ex ex' := match (ex' - ex)%Z with | Zneg d => (shift_pos d mx, ex') @@ -1049,7 +1130,8 @@ Theorem binary_round_correct : let x := F2R (Float radix2 (cond_Zopp sx (Zpos mx)) ex) in if Rlt_bool (Rabs (round radix2 fexp (round_mode m) x)) (bpow radix2 emax) then FF2R radix2 z = round radix2 fexp (round_mode m) x /\ - is_finite_FF z = true + is_finite_FF z = true /\ + sign_FF z = sx else z = binary_overflow m sx. Proof. @@ -1084,22 +1166,35 @@ Theorem binary_normalize_correct : forall m mx ex szero, if Rlt_bool (Rabs (round radix2 fexp (round_mode m) (F2R (Float radix2 mx ex)))) (bpow radix2 emax) then B2R (binary_normalize m mx ex szero) = round radix2 fexp (round_mode m) (F2R (Float radix2 mx ex)) /\ - is_finite (binary_normalize m mx ex szero) = true + is_finite (binary_normalize m mx ex szero) = true /\ + Bsign (binary_normalize m mx ex szero) = + match Rcompare (F2R (Float radix2 mx ex)) 0 with + | Eq => szero + | Lt => true + | Gt => false + end else B2FF (binary_normalize m mx ex szero) = binary_overflow m (Rlt_bool (F2R (Float radix2 mx ex)) 0). Proof with auto with typeclass_instances. intros m mx ez szero. destruct mx as [|mz|mz] ; simpl. rewrite F2R_0, round_0, Rabs_R0, Rlt_bool_true... +split... split... +rewrite Rcompare_Eq... apply bpow_gt_0. (* . mz > 0 *) generalize (binary_round_correct m false mz ez). simpl. case Rlt_bool_spec. -intros _ (Vz, (Rz, Rz')). +intros _ (Vz, (Rz, (Rz', Rz''))). split. now rewrite B2R_FF2B. +split. now rewrite is_finite_FF2B. +rewrite Bsign_FF2B, Rz''. +rewrite Rcompare_Gt... +apply F2R_gt_0_compat. +simpl. zify; omega. intros Hz' (Vz, Rz). rewrite B2FF_FF2B, Rz. apply f_equal. @@ -1110,10 +1205,15 @@ now apply F2R_ge_0_compat. generalize (binary_round_correct m true mz ez). simpl. case Rlt_bool_spec. -intros _ (Vz, (Rz, Rz')). +intros _ (Vz, (Rz, (Rz', Rz''))). split. now rewrite B2R_FF2B. +split. now rewrite is_finite_FF2B. +rewrite Bsign_FF2B, Rz''. +rewrite Rcompare_Lt... +apply F2R_lt_0_compat. +simpl. zify; omega. intros Hz' (Vz, Rz). rewrite B2FF_FF2B, Rz. apply f_equal. @@ -1123,12 +1223,12 @@ now apply F2R_lt_0_compat. Qed. (** Addition *) -Definition Bplus m x y := +Definition Bplus plus_nan m x y := + let f pl := B754_nan (fst pl) (snd pl) in match x, y with - | B754_nan, _ => x - | _, B754_nan => y + | B754_nan _ _, _ | _, B754_nan _ _ => f (plus_nan x y) | B754_infinity sx, B754_infinity sy => - if Bool.eqb sx sy then x else B754_nan + if Bool.eqb sx sy then x else f (plus_nan x y) | B754_infinity _, _ => x | _, B754_infinity _ => y | B754_zero sx, B754_zero sy => @@ -1143,28 +1243,47 @@ Definition Bplus m x y := end. Theorem Bplus_correct : - forall m x y, + forall plus_nan m x y, is_finite x = true -> is_finite y = true -> if Rlt_bool (Rabs (round radix2 fexp (round_mode m) (B2R x + B2R y))) (bpow radix2 emax) then - B2R (Bplus m x y) = round radix2 fexp (round_mode m) (B2R x + B2R y) /\ - is_finite (Bplus m x y) = true + B2R (Bplus plus_nan m x y) = round radix2 fexp (round_mode m) (B2R x + B2R y) /\ + is_finite (Bplus plus_nan m x y) = true /\ + Bsign (Bplus plus_nan m x y) = + match Rcompare (B2R x + B2R y) 0 with + | Eq => match m with mode_DN => orb (Bsign x) (Bsign y) + | _ => andb (Bsign x) (Bsign y) end + | Lt => true + | Gt => false + end else - (B2FF (Bplus m x y) = binary_overflow m (Bsign x) /\ Bsign x = Bsign y). + (B2FF (Bplus plus_nan m x y) = binary_overflow m (Bsign x) /\ Bsign x = Bsign y). Proof with auto with typeclass_instances. -intros m [sx|sx| |sx mx ex Hx] [sy|sy| |sy my ey Hy] Fx Fy ; try easy. +intros plus_nan m [sx|sx| |sx mx ex Hx] [sy|sy| |sy my ey Hy] Fx Fy ; try easy. (* *) rewrite Rplus_0_r, round_0, Rabs_R0, Rlt_bool_true... simpl. -case (Bool.eqb sx sy) ; try easy. -now case m. +rewrite Rcompare_Eq by auto. +destruct sx, sy; try easy; now case m. apply bpow_gt_0. (* *) rewrite Rplus_0_l, round_generic, Rlt_bool_true... +split... split... +simpl. unfold F2R. +erewrite <- Rmult_0_l, Rcompare_mult_r. +rewrite Rcompare_Z2R with (y:=0%Z). +destruct sy... +apply bpow_gt_0. apply abs_B2R_lt_emax. apply generic_format_B2R. (* *) rewrite Rplus_0_r, round_generic, Rlt_bool_true... +split... split... +simpl. unfold F2R. +erewrite <- Rmult_0_l, Rcompare_mult_r. +rewrite Rcompare_Z2R with (y:=0%Z). +destruct sx... +apply bpow_gt_0. apply abs_B2R_lt_emax. apply generic_format_B2R. (* *) @@ -1264,7 +1383,19 @@ now apply F2R_le_0_compat. (* . *) generalize (binary_normalize_correct m mz ez szero). case Rlt_bool_spec. -easy. +split; try easy. split; try easy. +destruct (Rcompare_spec (F2R (beta:=radix2) {| Fnum := mz; Fexp := ez |}) 0); try easy. +rewrite H1 in Hp. +apply Rplus_opp_r_uniq in Hp. +rewrite <- F2R_Zopp in Hp. +eapply canonic_unicity in Hp. +inversion Hp. destruct sy, sx, m; try discriminate H3; easy. +apply canonic_canonic_mantissa. +apply Bool.andb_true_iff in Hy. easy. +replace (-cond_Zopp sx (Z.pos mx))%Z with (cond_Zopp (negb sx) (Z.pos mx)) + by (destruct sx; auto). +apply canonic_canonic_mantissa. +apply Bool.andb_true_iff in Hx. easy. intros Hz' Vz. specialize (Sz Hz'). split. @@ -1273,26 +1404,32 @@ now apply f_equal. apply Sz. Qed. -Definition Bminus m x y := Bplus m x (Bopp y). +Definition Bminus minus_nan m x y := Bplus minus_nan m x (Bopp pair y). Theorem Bminus_correct : - forall m x y, + forall minus_nan m x y, is_finite x = true -> is_finite y = true -> if Rlt_bool (Rabs (round radix2 fexp (round_mode m) (B2R x - B2R y))) (bpow radix2 emax) then - B2R (Bminus m x y) = round radix2 fexp (round_mode m) (B2R x - B2R y) /\ - is_finite (Bminus m x y) = true + B2R (Bminus minus_nan m x y) = round radix2 fexp (round_mode m) (B2R x - B2R y) /\ + is_finite (Bminus minus_nan m x y) = true /\ + Bsign (Bminus minus_nan m x y) = + match Rcompare (B2R x - B2R y) 0 with + | Eq => match m with mode_DN => orb (Bsign x) (negb (Bsign y)) + | _ => andb (Bsign x) (negb (Bsign y)) end + | Lt => true + | Gt => false + end else - (B2FF (Bminus m x y) = binary_overflow m (Bsign x) /\ Bsign x = negb (Bsign y)). + (B2FF (Bminus minus_nan m x y) = binary_overflow m (Bsign x) /\ Bsign x = negb (Bsign y)). Proof with auto with typeclass_instances. -intros m x y Fx Fy. -replace (negb (Bsign y)) with (Bsign (Bopp y)). +intros m minus_nan x y Fx Fy. +replace (negb (Bsign y)) with (Bsign (Bopp pair y)). unfold Rminus. -rewrite <- B2R_Bopp. +erewrite <- B2R_Bopp. apply Bplus_correct. exact Fx. -now rewrite is_finite_Bopp. -now destruct y as [ | | | ]. +rewrite is_finite_Bopp. auto. now destruct y as [ | | | ]. Qed. (** Division *) @@ -1316,12 +1453,12 @@ Lemma Bdiv_correct_aux : let '(mz, ez, lz) := Fdiv_core_binary (Zpos mx) ex (Zpos my) ey in match mz with | Zpos mz => binary_round_aux m (xorb sx sy) mz ez lz - | _ => F754_nan (* dummy *) + | _ => F754_nan false xH (* dummy *) end in valid_binary z = true /\ if Rlt_bool (Rabs (round radix2 fexp (round_mode m) (x / y))) (bpow radix2 emax) then FF2R radix2 z = round radix2 fexp (round_mode m) (x / y) /\ - is_finite_FF z = true + is_finite_FF z = true /\ sign_FF z = xorb sx sy else z = binary_overflow m (xorb sx sy). Proof. @@ -1406,45 +1543,49 @@ apply Rinv_0_lt_compat. now apply F2R_gt_0_compat. Qed. -Definition Bdiv m x y := +Definition Bdiv div_nan m x y := + let f pl := B754_nan (fst pl) (snd pl) in match x, y with - | B754_nan, _ => x - | _, B754_nan => y - | B754_infinity sx, B754_infinity sy => B754_nan + | B754_nan _ _, _ | _, B754_nan _ _ => f (div_nan x y) + | B754_infinity sx, B754_infinity sy => f (div_nan x y) | B754_infinity sx, B754_finite sy _ _ _ => B754_infinity (xorb sx sy) | B754_finite sx _ _ _, B754_infinity sy => B754_zero (xorb sx sy) | B754_infinity sx, B754_zero sy => B754_infinity (xorb sx sy) | B754_zero sx, B754_infinity sy => B754_zero (xorb sx sy) | B754_finite sx _ _ _, B754_zero sy => B754_infinity (xorb sx sy) | B754_zero sx, B754_finite sy _ _ _ => B754_zero (xorb sx sy) - | B754_zero sx, B754_zero sy => B754_nan + | B754_zero sx, B754_zero sy => f (div_nan x y) | B754_finite sx mx ex _, B754_finite sy my ey _ => FF2B _ (proj1 (Bdiv_correct_aux m sx mx ex sy my ey)) end. Theorem Bdiv_correct : - forall m x y, + forall div_nan m x y, B2R y <> R0 -> if Rlt_bool (Rabs (round radix2 fexp (round_mode m) (B2R x / B2R y))) (bpow radix2 emax) then - B2R (Bdiv m x y) = round radix2 fexp (round_mode m) (B2R x / B2R y) /\ - is_finite (Bdiv m x y) = is_finite x + B2R (Bdiv div_nan m x y) = round radix2 fexp (round_mode m) (B2R x / B2R y) /\ + is_finite (Bdiv div_nan m x y) = is_finite x /\ + (is_nan (Bdiv div_nan m x y) = false -> + Bsign (Bdiv div_nan m x y) = xorb (Bsign x) (Bsign y)) else - B2FF (Bdiv m x y) = binary_overflow m (xorb (Bsign x) (Bsign y)). + B2FF (Bdiv div_nan m x y) = binary_overflow m (xorb (Bsign x) (Bsign y)). Proof. -intros m x [sy|sy| |sy my ey Hy] Zy ; try now elim Zy. +intros div_nan m x [sy|sy|sy ply|sy my ey Hy] Zy ; try now elim Zy. revert x. unfold Rdiv. -intros [sx|sx| |sx mx ex Hx] ; - try ( rewrite Rmult_0_l, round_0, Rabs_R0, Rlt_bool_true ; [ split ; apply refl_equal | apply bpow_gt_0 | auto with typeclass_instances ] ). +intros [sx|sx|sx plx|sx mx ex Hx] ; + try ( rewrite Rmult_0_l, round_0, Rabs_R0, Rlt_bool_true ; [ now repeat constructor | apply bpow_gt_0 | auto with typeclass_instances ] ). simpl. case Bdiv_correct_aux. intros H1. unfold Rdiv. case Rlt_bool. -intros (H2, H3). +intros (H2, (H3, H4)). split. now rewrite B2R_FF2B. +split. now rewrite is_finite_FF2B. +rewrite Bsign_FF2B. congruence. intros H2. now rewrite B2FF_FF2B. Qed. @@ -1473,11 +1614,11 @@ Lemma Bsqrt_correct_aux : let '(mz, ez, lz) := Fsqrt_core_binary (Zpos mx) ex in match mz with | Zpos mz => binary_round_aux m false mz ez lz - | _ => F754_nan (* dummy *) + | _ => F754_nan false xH (* dummy *) end in valid_binary z = true /\ FF2R radix2 z = round radix2 fexp (round_mode m) (sqrt x) /\ - is_finite_FF z = true. + is_finite_FF z = true /\ sign_FF z = false. Proof with auto with typeclass_instances. intros m mx ex Hx. replace (Fsqrt_core_binary (Zpos mx) ex) with (Fsqrt_core radix2 prec (Zpos mx) ex). @@ -1578,28 +1719,30 @@ now case mz. apply sqrt_ge_0. Qed. -Definition Bsqrt m x := +Definition Bsqrt sqrt_nan m x := + let f pl := B754_nan (fst pl) (snd pl) in match x with - | B754_nan => x + | B754_nan sx plx => f (sqrt_nan x) | B754_infinity false => x - | B754_infinity true => B754_nan - | B754_finite true _ _ _ => B754_nan + | B754_infinity true => f (sqrt_nan x) + | B754_finite true _ _ _ => f (sqrt_nan x) | B754_zero _ => x | B754_finite sx mx ex Hx => FF2B _ (proj1 (Bsqrt_correct_aux m mx ex Hx)) end. Theorem Bsqrt_correct : - forall m x, - B2R (Bsqrt m x) = round radix2 fexp (round_mode m) (sqrt (B2R x)) /\ - is_finite (Bsqrt m x) = match x with B754_zero _ => true | B754_finite false _ _ _ => true | _ => false end. + forall sqrt_nan m x, + B2R (Bsqrt sqrt_nan m x) = round radix2 fexp (round_mode m) (sqrt (B2R x)) /\ + is_finite (Bsqrt sqrt_nan m x) = match x with B754_zero _ => true | B754_finite false _ _ _ => true | _ => false end /\ + (is_nan (Bsqrt sqrt_nan m x) = false -> Bsign (Bsqrt sqrt_nan m x) = Bsign x). Proof. -intros m [sx|[|]| |sx mx ex Hx] ; try ( now simpl ; rewrite sqrt_0, round_0 ; auto with typeclass_instances ). +intros sqrt_nan m [sx|[|]| |sx mx ex Hx] ; try ( now simpl ; rewrite sqrt_0, round_0 ; intuition auto with typeclass_instances ). simpl. case Bsqrt_correct_aux. -intros H1 (H2, H3). +intros H1 (H2, (H3, H4)). case sx. -refine (conj _ (refl_equal false)). +refine (conj _ (conj (refl_equal false) _)). apply sym_eq. unfold sqrt. case Rcase_abs. @@ -1609,9 +1752,12 @@ auto with typeclass_instances. intros H. elim Rge_not_lt with (1 := H). now apply F2R_lt_0_compat. +easy. split. now rewrite B2R_FF2B. +split. now rewrite is_finite_FF2B. +intro. rewrite Bsign_FF2B. auto. Qed. End Binary. diff --git a/flocq/Appli/Fappli_IEEE_bits.v b/flocq/Appli/Fappli_IEEE_bits.v index 06ed21e..a41fba9 100644 --- a/flocq/Appli/Fappli_IEEE_bits.v +++ b/flocq/Appli/Fappli_IEEE_bits.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2011 Sylvie Boldo +Copyright (C) 2011-2013 Sylvie Boldo #
# -Copyright (C) 2011 Guillaume Melquiond +Copyright (C) 2011-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public @@ -172,7 +172,7 @@ Definition bits_of_binary_float (x : binary_float) := match x with | B754_zero sx => join_bits sx 0 0 | B754_infinity sx => join_bits sx 0 (Zpower 2 ew - 1) - | B754_nan => join_bits false (Zpower 2 mw - 1) (Zpower 2 ew - 1) + | B754_nan sx (exist plx _) => join_bits sx (Zpos plx) (Zpower 2 ew - 1) | B754_finite sx mx ex _ => if Zle_bool (Zpower 2 mw) (Zpos mx) then join_bits sx (Zpos mx - Zpower 2 mw) (ex - emin + 1) @@ -184,7 +184,7 @@ Definition split_bits_of_binary_float (x : binary_float) := match x with | B754_zero sx => (sx, 0, 0)%Z | B754_infinity sx => (sx, 0, Zpower 2 ew - 1)%Z - | B754_nan => (false, Zpower 2 mw - 1, Zpower 2 ew - 1)%Z + | B754_nan sx (exist plx _) => (sx, Zpos plx, Zpower 2 ew - 1)%Z | B754_finite sx mx ex _ => if Zle_bool (Zpower 2 mw) (Zpos mx) then (sx, Zpos mx - Zpower 2 mw, ex - emin + 1)%Z @@ -196,8 +196,16 @@ Theorem split_bits_of_binary_float_correct : forall x, split_bits (bits_of_binary_float x) = split_bits_of_binary_float x. Proof. -intros [sx|sx| |sx mx ex Hx] ; +intros [sx|sx|sx [plx Hplx]|sx mx ex Hx] ; try ( simpl ; apply split_join_bits ; split ; try apply Zle_refl ; try apply Zlt_pred ; trivial ; omega ). +simpl. apply split_join_bits; split; try (zify; omega). +destruct (digits2_Pnat_correct plx). +rewrite Zpower_nat_Z in H0. +eapply Zlt_le_trans. apply H0. +change 2%Z with (radix_val radix2). apply Zpower_le. +rewrite Z.ltb_lt in Hplx. +unfold prec in *. zify; omega. +(* *) unfold bits_of_binary_float, split_bits_of_binary_float. assert (Hf: (emin <= ex /\ Zdigits radix2 (Zpos mx) <= prec)%Z). destruct (andb_prop _ _ Hx) as (Hx', _). @@ -246,14 +254,18 @@ Definition binary_float_of_bits_aux x := match mx with | Z0 => F754_zero sx | Zpos px => F754_finite sx px emin - | Zneg _ => F754_nan (* dummy *) + | Zneg _ => F754_nan false xH (* dummy *) end else if Zeq_bool ex (Zpower 2 ew - 1) then - if Zeq_bool mx 0 then F754_infinity sx else F754_nan + match mx with + | Z0 => F754_infinity sx + | Zpos plx => F754_nan sx plx + | Zneg _ => F754_nan false xH (* dummy *) + end else match (mx + Zpower 2 mw)%Z with | Zpos px => F754_finite sx px (ex + emin - 1) - | _ => F754_nan (* dummy *) + | _ => F754_nan false xH (* dummy *) end. Lemma binary_float_of_bits_aux_correct : @@ -292,9 +304,20 @@ cut (0 < emax)%Z. clear -H Hew ; omega. apply (Zpower_gt_0 radix2). clear -Hew ; omega. apply bpow_gt_0. +simpl. intros. rewrite Z.ltb_lt. unfold prec. zify; omega. case Zeq_bool_spec ; intros He2. -now case Zeq_bool. +case_eq (x mod 2 ^ mw)%Z; try easy. +(* nan *) +intros plx Eqplx. apply Z.ltb_lt. +rewrite Z_of_nat_S_digits2_Pnat. +assert (forall a b, a <= b -> a < b+1)%Z by (intros; omega). apply H. clear H. +apply Zdigits_le_Zpower. simpl. +rewrite <- Eqplx. edestruct Z_mod_lt; eauto. +change 2%Z with (radix_val radix2). +apply Z.lt_gt, Zpower_gt_0. omega. +simpl. intros. rewrite Z.ltb_lt. unfold prec. zify; omega. case_eq (x mod 2^mw + 2^mw)%Z ; try easy. +simpl. intros. rewrite Z.ltb_lt. unfold prec. zify; omega. (* normal *) intros px Hm. assert (prec = Zdigits radix2 (Zpos px)). @@ -365,6 +388,7 @@ apply Zlt_gt. apply (Zpower_gt_0 radix2). now apply Zlt_le_weak. apply bpow_gt_0. +simpl. intros. rewrite Z.ltb_lt. unfold prec. zify; omega. Qed. Definition binary_float_of_bits x := @@ -380,7 +404,7 @@ unfold binary_float_of_bits. rewrite B2FF_FF2B. unfold binary_float_of_bits_aux. rewrite split_bits_of_binary_float_correct. -destruct x as [sx|sx| |sx mx ex Bx]. +destruct x as [sx|sx|sx [plx Hplx]|sx mx ex Bx]. apply refl_equal. (* *) simpl. @@ -391,12 +415,7 @@ now apply (Zpower_gt_1 radix2). (* *) simpl. rewrite Zeq_bool_false. -rewrite Zeq_bool_true. -rewrite Zeq_bool_false. -apply refl_equal. -cut (1 < 2^mw)%Z. clear ; omega. -now apply (Zpower_gt_1 radix2). -apply refl_equal. +rewrite Zeq_bool_true; auto. cut (1 < 2^ew)%Z. clear ; omega. now apply (Zpower_gt_1 radix2). (* *) @@ -442,7 +461,6 @@ Qed. Theorem bits_of_binary_float_of_bits : forall x, (0 <= x < 2^(mw+ew+1))%Z -> - binary_float_of_bits x <> B754_nan prec emax -> bits_of_binary_float (binary_float_of_bits x) = x. Proof. intros x Hx. @@ -462,28 +480,28 @@ now apply Zlt_gt. case Zeq_bool_spec ; intros He1. (* subnormal *) case_eq mx. -intros Hm Jx _ _. +intros Hm Jx _. now rewrite He1 in Jx. -intros px Hm Jx _ _. +intros px Hm Jx _. rewrite Zle_bool_false. now rewrite <- He1. now rewrite <- Hm. -intros px Hm _ _ _. +intros px Hm _ _. apply False_ind. apply Zle_not_lt with (1 := proj1 Bm). now rewrite Hm. case Zeq_bool_spec ; intros He2. (* infinity/nan *) -case Zeq_bool_spec ; intros Hm. -now rewrite Hm, He2. -intros _ Cx Nx. -now elim Nx. +case_eq mx; intros Hm. +now rewrite He2. +now rewrite He2. +intros. zify; omega. (* normal *) case_eq (mx + 2 ^ mw)%Z. intros Hm. apply False_ind. clear -Bm Hm ; omega. -intros p Hm Jx Cx _. +intros p Hm Jx Cx. rewrite <- Hm. rewrite Zle_bool_true. now ring_simplify (mx + 2^mw - 2^mw)%Z (ex + emin - 1 - emin + 1)%Z. diff --git a/flocq/Appli/Fappli_rnd_odd.v b/flocq/Appli/Fappli_rnd_odd.v new file mode 100644 index 0000000..b4a2c52 --- /dev/null +++ b/flocq/Appli/Fappli_rnd_odd.v @@ -0,0 +1,979 @@ +(** +This file is part of the Flocq formalization of floating-point +arithmetic in Coq: http://flocq.gforge.inria.fr/ + +Copyright (C) 2010-2013 Sylvie Boldo +#
# +Copyright (C) 2010-2013 Guillaume Melquiond + +This library is free software; you can redistribute it and/or +modify it under the terms of the GNU Lesser General Public +License as published by the Free Software Foundation; either +version 3 of the License, or (at your option) any later version. + +This library is distributed in the hope that it will be useful, +but WITHOUT ANY WARRANTY; without even the implied warranty of +MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +COPYING file for more details. +*) + +(** * Rounding to odd and its properties, including the equivalence + between rnd_NE and double rounding with rnd_odd and then rnd_NE *) + + +Require Import Fcore. +Require Import Fcalc_ops. + +Definition Zrnd_odd x := match Req_EM_T x (Z2R (Zfloor x)) with + | left _ => Zfloor x + | right _ => match (Zeven (Zfloor x)) with + | true => Zceil x + | false => Zfloor x + end + end. + + + +Global Instance valid_rnd_odd : Valid_rnd Zrnd_odd. +Proof. +split. +(* . *) +intros x y Hxy. +assert (Zfloor x <= Zrnd_odd y)%Z. +(* .. *) +apply Zle_trans with (Zfloor y). +now apply Zfloor_le. +unfold Zrnd_odd; destruct (Req_EM_T y (Z2R (Zfloor y))). +now apply Zle_refl. +case (Zeven (Zfloor y)). +apply le_Z2R. +apply Rle_trans with y. +apply Zfloor_lb. +apply Zceil_ub. +now apply Zle_refl. +unfold Zrnd_odd at 1. +(* . *) +destruct (Req_EM_T x (Z2R (Zfloor x))) as [Hx|Hx]. +(* .. *) +apply H. +(* .. *) +case_eq (Zeven (Zfloor x)); intros Hx2. +2: apply H. +unfold Zrnd_odd; destruct (Req_EM_T y (Z2R (Zfloor y))) as [Hy|Hy]. +apply Zceil_glb. +now rewrite <- Hy. +case_eq (Zeven (Zfloor y)); intros Hy2. +now apply Zceil_le. +apply Zceil_glb. +assert (H0:(Zfloor x <= Zfloor y)%Z) by now apply Zfloor_le. +case (Zle_lt_or_eq _ _ H0); intros H1. +apply Rle_trans with (1:=Zceil_ub _). +rewrite Zceil_floor_neq. +apply Z2R_le; omega. +now apply sym_not_eq. +contradict Hy2. +rewrite <- H1, Hx2; discriminate. +(* . *) +intros n; unfold Zrnd_odd. +rewrite Zfloor_Z2R, Zceil_Z2R. +destruct (Req_EM_T (Z2R n) (Z2R n)); trivial. +case (Zeven n); trivial. +Qed. + + + +Lemma Zrnd_odd_Zodd: forall x, x <> (Z2R (Zfloor x)) -> + (Zeven (Zrnd_odd x)) = false. +Proof. +intros x Hx; unfold Zrnd_odd. +destruct (Req_EM_T x (Z2R (Zfloor x))) as [H|H]. +now contradict H. +case_eq (Zeven (Zfloor x)). +(* difficult case *) +intros H'. +rewrite Zceil_floor_neq. +rewrite Zeven_plus, H'. +reflexivity. +now apply sym_not_eq. +trivial. +Qed. + + + + +Section Fcore_rnd_odd. + +Variable beta : radix. + +Notation bpow e := (bpow beta e). + +Variable fexp : Z -> Z. + +Context { valid_exp : Valid_exp fexp }. +Context { exists_NE_ : Exists_NE beta fexp }. + +Notation format := (generic_format beta fexp). +Notation canonic := (canonic beta fexp). +Notation cexp := (canonic_exp beta fexp). + + +Definition Rnd_odd_pt (x f : R) := + format f /\ ((f = x)%R \/ + ((Rnd_DN_pt format x f \/ Rnd_UP_pt format x f) /\ + exists g : float beta, f = F2R g /\ canonic g /\ Zeven (Fnum g) = false)). + +Definition Rnd_odd (rnd : R -> R) := + forall x : R, Rnd_odd_pt x (rnd x). + + +Theorem Rnd_odd_pt_sym : forall x f : R, + Rnd_odd_pt (-x) (-f) -> Rnd_odd_pt x f. +Proof with auto with typeclass_instances. +intros x f (H1,H2). +split. +replace f with (-(-f))%R by ring. +now apply generic_format_opp. +destruct H2. +left. +replace f with (-(-f))%R by ring. +rewrite H; ring. +right. +destruct H as (H2,(g,(Hg1,(Hg2,Hg3)))). +split. +destruct H2. +right. +replace f with (-(-f))%R by ring. +replace x with (-(-x))%R by ring. +apply Rnd_DN_UP_pt_sym... +apply generic_format_opp. +left. +replace f with (-(-f))%R by ring. +replace x with (-(-x))%R by ring. +apply Rnd_UP_DN_pt_sym... +apply generic_format_opp. +exists (Float beta (-Fnum g) (Fexp g)). +split. +rewrite F2R_Zopp. +replace f with (-(-f))%R by ring. +rewrite Hg1; reflexivity. +split. +now apply canonic_opp. +simpl. +now rewrite Zeven_opp. +Qed. + + +Theorem round_odd_opp : + forall x, + (round beta fexp Zrnd_odd (-x) = (- round beta fexp Zrnd_odd x))%R. +Proof. +intros x; unfold round. +rewrite <- F2R_Zopp. +unfold F2R; simpl. +apply f_equal2; apply f_equal. +rewrite scaled_mantissa_opp. +generalize (scaled_mantissa beta fexp x); intros r. +unfold Zrnd_odd. +case (Req_EM_T (- r) (Z2R (Zfloor (- r)))). +case (Req_EM_T r (Z2R (Zfloor r))). +intros Y1 Y2. +apply eq_Z2R. +now rewrite Z2R_opp, <- Y1, <-Y2. +intros Y1 Y2. +absurd (r=Z2R (Zfloor r)); trivial. +pattern r at 2; replace r with (-(-r))%R by ring. +rewrite Y2, <- Z2R_opp. +rewrite Zfloor_Z2R. +rewrite Z2R_opp, <- Y2. +ring. +case (Req_EM_T r (Z2R (Zfloor r))). +intros Y1 Y2. +absurd (-r=Z2R (Zfloor (-r)))%R; trivial. +pattern r at 2; rewrite Y1. +rewrite <- Z2R_opp, Zfloor_Z2R. +now rewrite Z2R_opp, <- Y1. +intros Y1 Y2. +unfold Zceil; rewrite Ropp_involutive. +replace (Zeven (Zfloor (- r))) with (negb (Zeven (Zfloor r))). +case (Zeven (Zfloor r)); simpl; ring. +apply trans_eq with (Zeven (Zceil r)). +rewrite Zceil_floor_neq. +rewrite Zeven_plus. +simpl; reflexivity. +now apply sym_not_eq. +rewrite <- (Zeven_opp (Zfloor (- r))). +reflexivity. +apply canonic_exp_opp. +Qed. + + + +Theorem round_odd_pt : + forall x, + Rnd_odd_pt x (round beta fexp Zrnd_odd x). +Proof with auto with typeclass_instances. +cut (forall x : R, (0 < x)%R -> Rnd_odd_pt x (round beta fexp Zrnd_odd x)). +intros H x; case (Rle_or_lt 0 x). +intros H1; destruct H1. +now apply H. +rewrite <- H0. +rewrite round_0... +split. +apply generic_format_0. +now left. +intros Hx; apply Rnd_odd_pt_sym. +rewrite <- round_odd_opp. +apply H. +auto with real. +(* *) +intros x Hxp. +generalize (generic_format_round beta fexp Zrnd_odd x). +set (o:=round beta fexp Zrnd_odd x). +intros Ho. +split. +assumption. +(* *) +case (Req_dec o x); intros Hx. +now left. +right. +assert (o=round beta fexp Zfloor x \/ o=round beta fexp Zceil x). +unfold o, round, F2R;simpl. +case (Zrnd_DN_or_UP Zrnd_odd (scaled_mantissa beta fexp x))... +intros H; rewrite H; now left. +intros H; rewrite H; now right. +split. +destruct H; rewrite H. +left; apply round_DN_pt... +right; apply round_UP_pt... +(* *) +unfold o, Zrnd_odd, round. +case (Req_EM_T (scaled_mantissa beta fexp x) + (Z2R (Zfloor (scaled_mantissa beta fexp x)))). +intros T. +absurd (o=x); trivial. +apply round_generic... +unfold generic_format, F2R; simpl. +rewrite <- (scaled_mantissa_mult_bpow beta fexp) at 1. +apply f_equal2; trivial; rewrite T at 1. +apply f_equal, sym_eq, Ztrunc_floor. +apply Rmult_le_pos. +now left. +apply bpow_ge_0. +intros L. +case_eq (Zeven (Zfloor (scaled_mantissa beta fexp x))). +(* . *) +generalize (generic_format_round beta fexp Zceil x). +unfold generic_format. +set (f:=round beta fexp Zceil x). +set (ef := canonic_exp beta fexp f). +set (mf := Ztrunc (scaled_mantissa beta fexp f)). +exists (Float beta mf ef). +unfold Fcore_generic_fmt.canonic. +rewrite <- H0. +repeat split; try assumption. +apply trans_eq with (negb (Zeven (Zfloor (scaled_mantissa beta fexp x)))). +2: rewrite H1; reflexivity. +apply trans_eq with (negb (Zeven (Fnum + (Float beta (Zfloor (scaled_mantissa beta fexp x)) (cexp x))))). +2: reflexivity. +case (Rle_lt_or_eq_dec 0 (round beta fexp Zfloor x)). +rewrite <- round_0 with beta fexp Zfloor... +apply round_le... +now left. +intros Y. +generalize (DN_UP_parity_generic beta fexp)... +unfold DN_UP_parity_prop. +intros T; apply T with x; clear T. +unfold generic_format. +rewrite <- (scaled_mantissa_mult_bpow beta fexp x) at 1. +unfold F2R; simpl. +apply Rmult_neq_compat_r. +apply Rgt_not_eq, bpow_gt_0. +rewrite Ztrunc_floor. +assumption. +apply Rmult_le_pos. +now left. +apply bpow_ge_0. +unfold Fcore_generic_fmt.canonic. +simpl. +apply sym_eq, canonic_exp_DN... +unfold Fcore_generic_fmt.canonic. +rewrite <- H0; reflexivity. +reflexivity. +apply trans_eq with (round beta fexp Ztrunc (round beta fexp Zceil x)). +reflexivity. +apply round_generic... +intros Y. +replace (Fnum {| Fnum := Zfloor (scaled_mantissa beta fexp x); Fexp := cexp x |}) + with (Fnum (Float beta 0 (fexp (ln_beta beta 0)))). +generalize (DN_UP_parity_generic beta fexp)... +unfold DN_UP_parity_prop. +intros T; apply T with x; clear T. +unfold generic_format. +rewrite <- (scaled_mantissa_mult_bpow beta fexp x) at 1. +unfold F2R; simpl. +apply Rmult_neq_compat_r. +apply Rgt_not_eq, bpow_gt_0. +rewrite Ztrunc_floor. +assumption. +apply Rmult_le_pos. +now left. +apply bpow_ge_0. +apply canonic_0. +unfold Fcore_generic_fmt.canonic. +rewrite <- H0; reflexivity. +rewrite <- Y; unfold F2R; simpl; ring. +apply trans_eq with (round beta fexp Ztrunc (round beta fexp Zceil x)). +reflexivity. +apply round_generic... +simpl. +apply eq_Z2R, Rmult_eq_reg_r with (bpow (cexp x)). +unfold round, F2R in Y; simpl in Y; rewrite <- Y. +simpl; ring. +apply Rgt_not_eq, bpow_gt_0. +(* . *) +intros Y. +case (Rle_lt_or_eq_dec 0 (round beta fexp Zfloor x)). +rewrite <- round_0 with beta fexp Zfloor... +apply round_le... +now left. +intros Hrx. +set (ef := canonic_exp beta fexp x). +set (mf := Zfloor (scaled_mantissa beta fexp x)). +exists (Float beta mf ef). +unfold Fcore_generic_fmt.canonic. +repeat split; try assumption. +simpl. +apply trans_eq with (cexp (round beta fexp Zfloor x )). +apply sym_eq, canonic_exp_DN... +reflexivity. +intros Hrx; contradict Y. +replace (Zfloor (scaled_mantissa beta fexp x)) with 0%Z. +simpl; discriminate. +apply eq_Z2R, Rmult_eq_reg_r with (bpow (cexp x)). +unfold round, F2R in Hrx; simpl in Hrx; rewrite <- Hrx. +simpl; ring. +apply Rgt_not_eq, bpow_gt_0. +Qed. + +End Fcore_rnd_odd. + +Section Odd_prop_aux. + +Variable beta : radix. +Hypothesis Even_beta: Zeven (radix_val beta)=true. + +Notation bpow e := (bpow beta e). + +Variable fexp : Z -> Z. +Variable fexpe : Z -> Z. + +Context { valid_exp : Valid_exp fexp }. +Context { exists_NE_ : Exists_NE beta fexp }. (* for underflow reason *) +Context { valid_expe : Valid_exp fexpe }. +Context { exists_NE_e : Exists_NE beta fexpe }. (* for defining rounding to odd *) + +Hypothesis fexpe_fexp: forall e, (fexpe e <= fexp e -2)%Z. + + +Lemma generic_format_fexpe_fexp: forall x, + generic_format beta fexp x -> generic_format beta fexpe x. +Proof. +intros x Hx. +apply generic_inclusion_ln_beta with fexp; trivial; intros Hx2. +generalize (fexpe_fexp (ln_beta beta x)). +omega. +Qed. + + + +Lemma exists_even_fexp_lt: forall (c:Z->Z), forall (x:R), + (exists f:float beta, F2R f = x /\ (c (ln_beta beta x) < Fexp f)%Z) -> + exists f:float beta, F2R f =x /\ canonic beta c f /\ Zeven (Fnum f) = true. +Proof with auto with typeclass_instances. +intros c x (g,(Hg1,Hg2)). +exists (Float beta + (Fnum g*Z.pow (radix_val beta) (Fexp g - c (ln_beta beta x))) + (c (ln_beta beta x))). +assert (F2R (Float beta + (Fnum g*Z.pow (radix_val beta) (Fexp g - c (ln_beta beta x))) + (c (ln_beta beta x))) = x). +unfold F2R; simpl. +rewrite Z2R_mult, Z2R_Zpower. +rewrite Rmult_assoc, <- bpow_plus. +rewrite <- Hg1; unfold F2R. +apply f_equal, f_equal. +ring. +omega. +split; trivial. +split. +unfold canonic, canonic_exp. +now rewrite H. +simpl. +rewrite Zeven_mult. +rewrite Zeven_Zpower. +rewrite Even_beta. +apply Bool.orb_true_intro. +now right. +omega. +Qed. + + +Variable choice:Z->bool. +Variable x:R. + + +Variable d u: float beta. +Hypothesis Hd: Rnd_DN_pt (generic_format beta fexp) x (F2R d). +Hypothesis Cd: canonic beta fexp d. +Hypothesis Hu: Rnd_UP_pt (generic_format beta fexp) x (F2R u). +Hypothesis Cu: canonic beta fexp u. + +Hypothesis xPos: (0 < x)%R. + + +Let m:= ((F2R d+F2R u)/2)%R. + + +Lemma d_eq: F2R d= round beta fexp Zfloor x. +Proof with auto with typeclass_instances. +apply Rnd_DN_pt_unicity with (generic_format beta fexp) x... +apply round_DN_pt... +Qed. + + +Lemma u_eq: F2R u= round beta fexp Zceil x. +Proof with auto with typeclass_instances. +apply Rnd_UP_pt_unicity with (generic_format beta fexp) x... +apply round_UP_pt... +Qed. + + +Lemma d_ge_0: (0 <= F2R d)%R. +Proof with auto with typeclass_instances. +rewrite d_eq; apply round_ge_generic... +apply generic_format_0. +now left. +Qed. + + + +Lemma ln_beta_d: (0< F2R d)%R -> + (ln_beta beta (F2R d) = ln_beta beta x :>Z). +Proof with auto with typeclass_instances. +intros Y. +rewrite d_eq; apply ln_beta_round_DN... +now rewrite <- d_eq. +Qed. + + +Lemma Fexp_d: (0 < F2R d)%R -> Fexp d =fexp (ln_beta beta x). +Proof with auto with typeclass_instances. +intros Y. +now rewrite Cd, <- ln_beta_d. +Qed. + + + +Lemma format_bpow_x: (0 < F2R d)%R + -> generic_format beta fexp (bpow (ln_beta beta x)). +Proof with auto with typeclass_instances. +intros Y. +apply generic_format_bpow. +apply valid_exp. +rewrite <- Fexp_d; trivial. +apply Zlt_le_trans with (ln_beta beta (F2R d))%Z. +rewrite Cd; apply ln_beta_generic_gt... +now apply Rgt_not_eq. +apply Hd. +apply ln_beta_le; trivial. +apply Hd. +Qed. + + +Lemma format_bpow_d: (0 < F2R d)%R -> + generic_format beta fexp (bpow (ln_beta beta (F2R d))). +Proof with auto with typeclass_instances. +intros Y; apply generic_format_bpow. +apply valid_exp. +apply ln_beta_generic_gt... +now apply Rgt_not_eq. +now apply generic_format_canonic. +Qed. + + +Lemma d_le_m: (F2R d <= m)%R. +apply Rmult_le_reg_l with 2%R. +auto with real. +apply Rplus_le_reg_l with (-F2R d)%R. +apply Rle_trans with (F2R d). +right; ring. +apply Rle_trans with (F2R u). +apply Rle_trans with x. +apply Hd. +apply Hu. +right; unfold m; field. +Qed. + +Lemma m_le_u: (m <= F2R u)%R. +apply Rmult_le_reg_l with 2%R. +auto with real. +apply Rplus_le_reg_l with (-F2R u)%R. +apply Rle_trans with (F2R d). +right; unfold m; field. +apply Rle_trans with (F2R u). +apply Rle_trans with x. +apply Hd. +apply Hu. +right; ring. +Qed. + +Lemma ln_beta_m: (0 < F2R d)%R -> (ln_beta beta m =ln_beta beta (F2R d) :>Z). +Proof with auto with typeclass_instances. +intros dPos; apply ln_beta_unique_pos. +split. +apply Rle_trans with (F2R d). +destruct (ln_beta beta (F2R d)) as (e,He). +simpl. +rewrite Rabs_right in He. +apply He. +now apply Rgt_not_eq. +apply Rle_ge; now left. +apply d_le_m. +case m_le_u; intros H. +apply Rlt_le_trans with (1:=H). +rewrite u_eq. +apply round_le_generic... +apply generic_format_bpow. +apply valid_exp. +apply ln_beta_generic_gt... +now apply Rgt_not_eq. +now apply generic_format_canonic. +case (Rle_or_lt x (bpow (ln_beta beta (F2R d)))); trivial; intros Z. +absurd ((bpow (ln_beta beta (F2R d)) <= (F2R d)))%R. +apply Rlt_not_le. +destruct (ln_beta beta (F2R d)) as (e,He). +simpl in *; rewrite Rabs_right in He. +apply He. +now apply Rgt_not_eq. +apply Rle_ge; now left. +apply Rle_trans with (round beta fexp Zfloor x). +2: right; apply sym_eq, d_eq. +apply round_ge_generic... +apply generic_format_bpow. +apply valid_exp. +apply ln_beta_generic_gt... +now apply Rgt_not_eq. +now apply generic_format_canonic. +now left. +replace m with (F2R d). +destruct (ln_beta beta (F2R d)) as (e,He). +simpl in *; rewrite Rabs_right in He. +apply He. +now apply Rgt_not_eq. +apply Rle_ge; now left. +assert (F2R d = F2R u). +apply Rmult_eq_reg_l with (/2)%R. +apply Rplus_eq_reg_l with (/2*F2R u)%R. +apply trans_eq with m. +unfold m, Rdiv; ring. +rewrite H; field. +auto with real. +apply Rgt_not_eq, Rlt_gt; auto with real. +unfold m; rewrite <- H0; field. +Qed. + + +Lemma ln_beta_m_0: (0 = F2R d)%R + -> (ln_beta beta m =ln_beta beta (F2R u)-1:>Z)%Z. +Proof with auto with typeclass_instances. +intros Y. +apply ln_beta_unique_pos. +unfold m; rewrite <- Y, Rplus_0_l. +rewrite u_eq. +destruct (ln_beta beta x) as (e,He). +rewrite Rabs_right in He. +rewrite round_UP_small_pos with (ex:=e). +rewrite ln_beta_bpow. +ring_simplify (fexp e + 1 - 1)%Z. +split. +unfold Zminus; rewrite bpow_plus. +unfold Rdiv; apply Rmult_le_compat_l. +apply bpow_ge_0. +simpl; unfold Z.pow_pos; simpl. +rewrite Zmult_1_r; apply Rinv_le. +auto with real. +apply (Z2R_le 2). +specialize (radix_gt_1 beta). +omega. +apply Rlt_le_trans with (bpow (fexp e)*1)%R. +2: right; ring. +unfold Rdiv; apply Rmult_lt_compat_l. +apply bpow_gt_0. +rewrite <- Rinv_1 at 3. +apply Rinv_lt; auto with real. +now apply He, Rgt_not_eq. +apply exp_small_round_0_pos with beta (Zfloor) x... +now apply He, Rgt_not_eq. +now rewrite <- d_eq, Y. +now left. +Qed. + + + + + +Lemma u'_eq: (0 < F2R d)%R -> exists f:float beta, F2R f = F2R u /\ (Fexp f = Fexp d)%Z. +Proof with auto with typeclass_instances. +intros Y. +rewrite u_eq; unfold round. +eexists; repeat split. +simpl; now rewrite Fexp_d. +Qed. + + + + +Lemma m_eq: (0 < F2R d)%R -> exists f:float beta, + F2R f = m /\ (Fexp f = fexp (ln_beta beta x) -1)%Z. +Proof with auto with typeclass_instances. +intros Y. +specialize (Zeven_ex (radix_val beta)); rewrite Even_beta. +intros (b, Hb); rewrite Zplus_0_r in Hb. +destruct u'_eq as (u', (Hu'1,Hu'2)); trivial. +exists (Fmult beta (Float beta b (-1)) (Fplus beta d u'))%R. +split. +rewrite F2R_mult, F2R_plus, Hu'1. +unfold m; rewrite Rmult_comm. +unfold Rdiv; apply f_equal. +unfold F2R; simpl; unfold Z.pow_pos; simpl. +rewrite Zmult_1_r, Hb, Z2R_mult. +simpl; field. +apply Rgt_not_eq, Rmult_lt_reg_l with (Z2R 2). +simpl; auto with real. +rewrite Rmult_0_r, <-Z2R_mult, <-Hb. +apply radix_pos. +apply trans_eq with (-1+Fexp (Fplus beta d u'))%Z. +unfold Fmult. +destruct (Fplus beta d u'); reflexivity. +rewrite Zplus_comm; unfold Zminus; apply f_equal2. +2: reflexivity. +rewrite Fexp_Fplus. +rewrite Z.min_l. +now rewrite Fexp_d. +rewrite Hu'2; omega. +Qed. + +Lemma m_eq_0: (0 = F2R d)%R -> exists f:float beta, + F2R f = m /\ (Fexp f = fexp (ln_beta beta (F2R u)) -1)%Z. +Proof with auto with typeclass_instances. +intros Y. +specialize (Zeven_ex (radix_val beta)); rewrite Even_beta. +intros (b, Hb); rewrite Zplus_0_r in Hb. +exists (Fmult beta (Float beta b (-1)) u)%R. +split. +rewrite F2R_mult; unfold m; rewrite <- Y, Rplus_0_l. +rewrite Rmult_comm. +unfold Rdiv; apply f_equal. +unfold F2R; simpl; unfold Z.pow_pos; simpl. +rewrite Zmult_1_r, Hb, Z2R_mult. +simpl; field. +apply Rgt_not_eq, Rmult_lt_reg_l with (Z2R 2). +simpl; auto with real. +rewrite Rmult_0_r, <-Z2R_mult, <-Hb. +apply radix_pos. +apply trans_eq with (-1+Fexp u)%Z. +unfold Fmult. +destruct u; reflexivity. +rewrite Zplus_comm, Cu; unfold Zminus; now apply f_equal2. +Qed. + +Lemma fexp_m_eq_0: (0 = F2R d)%R -> + (fexp (ln_beta beta (F2R u)-1) < fexp (ln_beta beta (F2R u))+1)%Z. +Proof with auto with typeclass_instances. +intros Y. +assert ((fexp (ln_beta beta (F2R u) - 1) <= fexp (ln_beta beta (F2R u))))%Z. +2: omega. +destruct (ln_beta beta x) as (e,He). +rewrite Rabs_right in He. +2: now left. +assert (e <= fexp e)%Z. +apply exp_small_round_0_pos with beta (Zfloor) x... +now apply He, Rgt_not_eq. +now rewrite <- d_eq, Y. +rewrite u_eq, round_UP_small_pos with (ex:=e); trivial. +2: now apply He, Rgt_not_eq. +rewrite ln_beta_bpow. +ring_simplify (fexp e + 1 - 1)%Z. +replace (fexp (fexp e)) with (fexp e). +case exists_NE_; intros V. +contradict V; rewrite Even_beta; discriminate. +rewrite (proj2 (V e)); omega. +apply sym_eq, valid_exp; omega. +Qed. + +Lemma Fm: generic_format beta fexpe m. +case (d_ge_0); intros Y. +(* *) +destruct m_eq as (g,(Hg1,Hg2)); trivial. +apply generic_format_F2R' with g. +now apply sym_eq. +intros H; unfold canonic_exp; rewrite Hg2. +rewrite ln_beta_m; trivial. +rewrite <- Fexp_d; trivial. +rewrite Cd. +unfold canonic_exp. +generalize (fexpe_fexp (ln_beta beta (F2R d))). +omega. +(* *) +destruct m_eq_0 as (g,(Hg1,Hg2)); trivial. +apply generic_format_F2R' with g. +assumption. +intros H; unfold canonic_exp; rewrite Hg2. +rewrite ln_beta_m_0; try assumption. +apply Zle_trans with (1:=fexpe_fexp _). +assert (fexp (ln_beta beta (F2R u)-1) < fexp (ln_beta beta (F2R u))+1)%Z;[idtac|omega]. +now apply fexp_m_eq_0. +Qed. + + + +Lemma Zm: + exists g : float beta, F2R g = m /\ canonic beta fexpe g /\ Zeven (Fnum g) = true. +Proof with auto with typeclass_instances. +case (d_ge_0); intros Y. +(* *) +destruct m_eq as (g,(Hg1,Hg2)); trivial. +apply exists_even_fexp_lt. +exists g; split; trivial. +rewrite Hg2. +rewrite ln_beta_m; trivial. +rewrite <- Fexp_d; trivial. +rewrite Cd. +unfold canonic_exp. +generalize (fexpe_fexp (ln_beta beta (F2R d))). +omega. +(* *) +destruct m_eq_0 as (g,(Hg1,Hg2)); trivial. +apply exists_even_fexp_lt. +exists g; split; trivial. +rewrite Hg2. +rewrite ln_beta_m_0; trivial. +apply Zle_lt_trans with (1:=fexpe_fexp _). +assert (fexp (ln_beta beta (F2R u)-1) < fexp (ln_beta beta (F2R u))+1)%Z;[idtac|omega]. +now apply fexp_m_eq_0. +Qed. + + +Lemma DN_odd_d_aux: forall z, (F2R d<= z< F2R u)%R -> + Rnd_DN_pt (generic_format beta fexp) z (F2R d). +Proof with auto with typeclass_instances. +intros z Hz1. +replace (F2R d) with (round beta fexp Zfloor z). +apply round_DN_pt... +case (Rnd_DN_UP_pt_split _ _ _ _ Hd Hu (round beta fexp Zfloor z)). +apply generic_format_round... +intros Y; apply Rle_antisym; trivial. +apply round_DN_pt... +apply Hd. +apply Hz1. +intros Y; absurd (z < z)%R. +auto with real. +apply Rlt_le_trans with (1:=proj2 Hz1), Rle_trans with (1:=Y). +apply round_DN_pt... +Qed. + +Lemma UP_odd_d_aux: forall z, (F2R d< z <= F2R u)%R -> + Rnd_UP_pt (generic_format beta fexp) z (F2R u). +Proof with auto with typeclass_instances. +intros z Hz1. +replace (F2R u) with (round beta fexp Zceil z). +apply round_UP_pt... +case (Rnd_DN_UP_pt_split _ _ _ _ Hd Hu (round beta fexp Zceil z)). +apply generic_format_round... +intros Y; absurd (z < z)%R. +auto with real. +apply Rle_lt_trans with (2:=proj1 Hz1), Rle_trans with (2:=Y). +apply round_UP_pt... +intros Y; apply Rle_antisym; trivial. +apply round_UP_pt... +apply Hu. +apply Hz1. +Qed. + + +Theorem round_odd_prop_pos: + round beta fexp (Znearest choice) (round beta fexpe Zrnd_odd x) = + round beta fexp (Znearest choice) x. +Proof with auto with typeclass_instances. +set (o:=round beta fexpe Zrnd_odd x). +case (generic_format_EM beta fexp x); intros Hx. +replace o with x; trivial. +unfold o; apply sym_eq, round_generic... +now apply generic_format_fexpe_fexp. +assert (K1:(F2R d <= o)%R). +apply round_ge_generic... +apply generic_format_fexpe_fexp, Hd. +apply Hd. +assert (K2:(o <= F2R u)%R). +apply round_le_generic... +apply generic_format_fexpe_fexp, Hu. +apply Hu. +assert (P:(x <> m -> o=m -> (forall P:Prop, P))). +intros Y1 Y2. +assert (H:(Rnd_odd_pt beta fexpe x o)). +apply round_odd_pt... +destruct H as (_,H); destruct H. +absurd (x=m)%R; try trivial. +now rewrite <- Y2, H. +destruct H as (_,(k,(Hk1,(Hk2,Hk3)))). +destruct Zm as (k',(Hk'1,(Hk'2,Hk'3))). +absurd (true=false). +discriminate. +rewrite <- Hk3, <- Hk'3. +apply f_equal, f_equal. +apply canonic_unicity with fexpe... +now rewrite Hk'1, <- Y2. +assert (generic_format beta fexp o -> (forall P:Prop, P)). +intros Y. +assert (H:(Rnd_odd_pt beta fexpe x o)). +apply round_odd_pt... +destruct H as (_,H); destruct H. +absurd (generic_format beta fexp x); trivial. +now rewrite <- H. +destruct H as (_,(k,(Hk1,(Hk2,Hk3)))). +destruct (exists_even_fexp_lt fexpe o) as (k',(Hk'1,(Hk'2,Hk'3))). +eexists; split. +apply sym_eq, Y. +simpl; unfold canonic_exp. +apply Zle_lt_trans with (1:=fexpe_fexp _). +omega. +absurd (true=false). +discriminate. +rewrite <- Hk3, <- Hk'3. +apply f_equal, f_equal. +apply canonic_unicity with fexpe... +now rewrite Hk'1, <- Hk1. +case K1; clear K1; intros K1. +2: apply H; rewrite <- K1; apply Hd. +case K2; clear K2; intros K2. +2: apply H; rewrite K2; apply Hu. +case (Rle_or_lt x m); intros Y;[destruct Y|idtac]. +(* . *) +apply trans_eq with (F2R d). +apply round_N_DN_betw with (F2R u)... +apply DN_odd_d_aux; split; try left; assumption. +apply UP_odd_d_aux; split; try left; assumption. +split. +apply round_ge_generic... +apply generic_format_fexpe_fexp, Hd. +apply Hd. +assert (o <= (F2R d + F2R u) / 2)%R. +apply round_le_generic... +apply Fm. +now left. +destruct H1; trivial. +apply P. +now apply Rlt_not_eq. +trivial. +apply sym_eq, round_N_DN_betw with (F2R u)... +split. +apply Hd. +exact H0. +(* . *) +replace o with x. +reflexivity. +apply sym_eq, round_generic... +rewrite H0; apply Fm. +(* . *) +apply trans_eq with (F2R u). +apply round_N_UP_betw with (F2R d)... +apply DN_odd_d_aux; split; try left; assumption. +apply UP_odd_d_aux; split; try left; assumption. +split. +assert ((F2R d + F2R u) / 2 <= o)%R. +apply round_ge_generic... +apply Fm. +now left. +destruct H0; trivial. +apply P. +now apply Rgt_not_eq. +rewrite <- H0; trivial. +apply round_le_generic... +apply generic_format_fexpe_fexp, Hu. +apply Hu. +apply sym_eq, round_N_UP_betw with (F2R d)... +split. +exact Y. +apply Hu. +Qed. + + +End Odd_prop_aux. + +Section Odd_prop. + +Variable beta : radix. +Hypothesis Even_beta: Zeven (radix_val beta)=true. + +Variable fexp : Z -> Z. +Variable fexpe : Z -> Z. +Variable choice:Z->bool. + +Context { valid_exp : Valid_exp fexp }. +Context { exists_NE_ : Exists_NE beta fexp }. (* for underflow reason *) +Context { valid_expe : Valid_exp fexpe }. +Context { exists_NE_e : Exists_NE beta fexpe }. (* for defining rounding to odd *) + +Hypothesis fexpe_fexp: forall e, (fexpe e <= fexp e -2)%Z. + + +Theorem canonizer: forall f, generic_format beta fexp f + -> exists g : float beta, f = F2R g /\ canonic beta fexp g. +Proof with auto with typeclass_instances. +intros f Hf. +exists (Float beta (Ztrunc (scaled_mantissa beta fexp f)) (canonic_exp beta fexp f)). +assert (L:(f = F2R (Float beta (Ztrunc (scaled_mantissa beta fexp f)) (canonic_exp beta fexp f)))). +apply trans_eq with (round beta fexp Ztrunc f). +apply sym_eq, round_generic... +reflexivity. +split; trivial. +unfold canonic; rewrite <- L. +reflexivity. +Qed. + + + + +Theorem round_odd_prop: forall x, + round beta fexp (Znearest choice) (round beta fexpe Zrnd_odd x) = + round beta fexp (Znearest choice) x. +Proof with auto with typeclass_instances. +intros x. +case (total_order_T x 0); intros H; [case H; clear H; intros H | idtac]. +rewrite <- (Ropp_involutive x). +rewrite round_odd_opp. +rewrite 2!round_N_opp. +apply f_equal. +destruct (canonizer (round beta fexp Zfloor (-x))) as (d,(Hd1,Hd2)). +apply generic_format_round... +destruct (canonizer (round beta fexp Zceil (-x))) as (u,(Hu1,Hu2)). +apply generic_format_round... +apply round_odd_prop_pos with d u... +rewrite <- Hd1; apply round_DN_pt... +rewrite <- Hu1; apply round_UP_pt... +auto with real. +(* . *) +rewrite H; repeat rewrite round_0... +(* . *) +destruct (canonizer (round beta fexp Zfloor x)) as (d,(Hd1,Hd2)). +apply generic_format_round... +destruct (canonizer (round beta fexp Zceil x)) as (u,(Hu1,Hu2)). +apply generic_format_round... +apply round_odd_prop_pos with d u... +rewrite <- Hd1; apply round_DN_pt... +rewrite <- Hu1; apply round_UP_pt... +Qed. + + +End Odd_prop. diff --git a/flocq/Calc/Fcalc_bracket.v b/flocq/Calc/Fcalc_bracket.v index dd4bd97..90a8588 100644 --- a/flocq/Calc/Fcalc_bracket.v +++ b/flocq/Calc/Fcalc_bracket.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2010-2011 Sylvie Boldo +Copyright (C) 2010-2013 Sylvie Boldo #
# -Copyright (C) 2010-2011 Guillaume Melquiond +Copyright (C) 2010-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public diff --git a/flocq/Calc/Fcalc_digits.v b/flocq/Calc/Fcalc_digits.v index 6210bac..4f76cc2 100644 --- a/flocq/Calc/Fcalc_digits.v +++ b/flocq/Calc/Fcalc_digits.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2010-2011 Sylvie Boldo +Copyright (C) 2010-2013 Sylvie Boldo #
# -Copyright (C) 2010-2011 Guillaume Melquiond +Copyright (C) 2010-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public diff --git a/flocq/Calc/Fcalc_div.v b/flocq/Calc/Fcalc_div.v index 6594e55..c8f1f9f 100644 --- a/flocq/Calc/Fcalc_div.v +++ b/flocq/Calc/Fcalc_div.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2010-2011 Sylvie Boldo +Copyright (C) 2010-2013 Sylvie Boldo #
# -Copyright (C) 2010-2011 Guillaume Melquiond +Copyright (C) 2010-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public diff --git a/flocq/Calc/Fcalc_ops.v b/flocq/Calc/Fcalc_ops.v index 15bb211..7ece683 100644 --- a/flocq/Calc/Fcalc_ops.v +++ b/flocq/Calc/Fcalc_ops.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2010-2011 Sylvie Boldo +Copyright (C) 2010-2013 Sylvie Boldo #
# -Copyright (C) 2010-2011 Guillaume Melquiond +Copyright (C) 2010-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public diff --git a/flocq/Calc/Fcalc_round.v b/flocq/Calc/Fcalc_round.v index 3d31aea..a1bcb84 100644 --- a/flocq/Calc/Fcalc_round.v +++ b/flocq/Calc/Fcalc_round.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2010-2011 Sylvie Boldo +Copyright (C) 2010-2013 Sylvie Boldo #
# -Copyright (C) 2010-2011 Guillaume Melquiond +Copyright (C) 2010-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public diff --git a/flocq/Calc/Fcalc_sqrt.v b/flocq/Calc/Fcalc_sqrt.v index e6fe74b..2ed3234 100644 --- a/flocq/Calc/Fcalc_sqrt.v +++ b/flocq/Calc/Fcalc_sqrt.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2010-2011 Sylvie Boldo +Copyright (C) 2010-2013 Sylvie Boldo #
# -Copyright (C) 2010-2011 Guillaume Melquiond +Copyright (C) 2010-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public diff --git a/flocq/Core/Fcore.v b/flocq/Core/Fcore.v index 23ebb39..2a5a5f0 100644 --- a/flocq/Core/Fcore.v +++ b/flocq/Core/Fcore.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2010-2011 Sylvie Boldo +Copyright (C) 2010-2013 Sylvie Boldo #
# -Copyright (C) 2010-2011 Guillaume Melquiond +Copyright (C) 2010-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public diff --git a/flocq/Core/Fcore_FIX.v b/flocq/Core/Fcore_FIX.v index f185ddf..a3b8d4d 100644 --- a/flocq/Core/Fcore_FIX.v +++ b/flocq/Core/Fcore_FIX.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2010-2011 Sylvie Boldo +Copyright (C) 2010-2013 Sylvie Boldo #
# -Copyright (C) 2010-2011 Guillaume Melquiond +Copyright (C) 2010-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public diff --git a/flocq/Core/Fcore_FLT.v b/flocq/Core/Fcore_FLT.v index 4ad4797..c057b6c 100644 --- a/flocq/Core/Fcore_FLT.v +++ b/flocq/Core/Fcore_FLT.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2010-2011 Sylvie Boldo +Copyright (C) 2010-2013 Sylvie Boldo #
# -Copyright (C) 2010-2011 Guillaume Melquiond +Copyright (C) 2010-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public diff --git a/flocq/Core/Fcore_FLX.v b/flocq/Core/Fcore_FLX.v index 62ecb6f..800540f 100644 --- a/flocq/Core/Fcore_FLX.v +++ b/flocq/Core/Fcore_FLX.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2010-2011 Sylvie Boldo +Copyright (C) 2010-2013 Sylvie Boldo #
# -Copyright (C) 2010-2011 Guillaume Melquiond +Copyright (C) 2010-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public diff --git a/flocq/Core/Fcore_FTZ.v b/flocq/Core/Fcore_FTZ.v index 5356c11..5f3e533 100644 --- a/flocq/Core/Fcore_FTZ.v +++ b/flocq/Core/Fcore_FTZ.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2010-2011 Sylvie Boldo +Copyright (C) 2010-2013 Sylvie Boldo #
# -Copyright (C) 2010-2011 Guillaume Melquiond +Copyright (C) 2010-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public diff --git a/flocq/Core/Fcore_Raux.v b/flocq/Core/Fcore_Raux.v index 748e36e..d811019 100644 --- a/flocq/Core/Fcore_Raux.v +++ b/flocq/Core/Fcore_Raux.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2010-2011 Sylvie Boldo +Copyright (C) 2010-2013 Sylvie Boldo #
# -Copyright (C) 2010-2011 Guillaume Melquiond +Copyright (C) 2010-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public @@ -121,6 +121,18 @@ rewrite <- 3!(Rmult_comm r). apply Rmult_minus_distr_l. Qed. +Lemma Rmult_neq_reg_r: forall r1 r2 r3:R, (r2 * r1 <> r3 * r1)%R -> r2 <> r3. +intros r1 r2 r3 H1 H2. +apply H1; rewrite H2; ring. +Qed. + +Lemma Rmult_neq_compat_r: forall r1 r2 r3:R, (r1 <> 0)%R -> (r2 <> r3)%R + -> (r2 *r1 <> r3*r1)%R. +intros r1 r2 r3 H H1 H2. +now apply H1, Rmult_eq_reg_r with r1. +Qed. + + Theorem Rmult_min_distr_r : forall r r1 r2 : R, (0 <= r)%R -> diff --git a/flocq/Core/Fcore_Zaux.v b/flocq/Core/Fcore_Zaux.v index af0d837..7ba28ca 100644 --- a/flocq/Core/Fcore_Zaux.v +++ b/flocq/Core/Fcore_Zaux.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2011 Sylvie Boldo +Copyright (C) 2011-2013 Sylvie Boldo #
# -Copyright (C) 2011 Guillaume Melquiond +Copyright (C) 2011-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public diff --git a/flocq/Core/Fcore_defs.v b/flocq/Core/Fcore_defs.v index fda3a85..ad8cf4f 100644 --- a/flocq/Core/Fcore_defs.v +++ b/flocq/Core/Fcore_defs.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2010-2011 Sylvie Boldo +Copyright (C) 2010-2013 Sylvie Boldo #
# -Copyright (C) 2010-2011 Guillaume Melquiond +Copyright (C) 2010-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public diff --git a/flocq/Core/Fcore_digits.v b/flocq/Core/Fcore_digits.v index 2ae076e..02c7a0e 100644 --- a/flocq/Core/Fcore_digits.v +++ b/flocq/Core/Fcore_digits.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2011 Sylvie Boldo +Copyright (C) 2011-2013 Sylvie Boldo #
# -Copyright (C) 2011 Guillaume Melquiond +Copyright (C) 2011-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public diff --git a/flocq/Core/Fcore_float_prop.v b/flocq/Core/Fcore_float_prop.v index 746f7a6..e1535bc 100644 --- a/flocq/Core/Fcore_float_prop.v +++ b/flocq/Core/Fcore_float_prop.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2010-2011 Sylvie Boldo +Copyright (C) 2010-2013 Sylvie Boldo #
# -Copyright (C) 2010-2011 Guillaume Melquiond +Copyright (C) 2010-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public diff --git a/flocq/Core/Fcore_generic_fmt.v b/flocq/Core/Fcore_generic_fmt.v index b1db47c..b04cf3d 100644 --- a/flocq/Core/Fcore_generic_fmt.v +++ b/flocq/Core/Fcore_generic_fmt.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2010-2011 Sylvie Boldo +Copyright (C) 2010-2013 Sylvie Boldo #
# -Copyright (C) 2010-2011 Guillaume Melquiond +Copyright (C) 2010-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public @@ -165,6 +165,22 @@ now rewrite Ztrunc_Z2R. now apply Zle_left. Qed. +Lemma generic_format_F2R': forall (x:R) (f:float beta), + F2R f = x -> ((x <> 0)%R -> + (canonic_exp x <= Fexp f)%Z) -> + generic_format x. +Proof. +intros x f H1 H2. +rewrite <- H1; destruct f as (m,e). +apply generic_format_F2R. +simpl in *; intros H3. +rewrite H1; apply H2. +intros Y; apply H3. +apply F2R_eq_0_reg with beta e. +now rewrite H1. +Qed. + + Theorem canonic_opp : forall m e, canonic (Float beta m e) -> @@ -175,6 +191,26 @@ unfold canonic. now rewrite F2R_Zopp, canonic_exp_opp. Qed. +Theorem canonic_abs : + forall m e, + canonic (Float beta m e) -> + canonic (Float beta (Zabs m) e). +Proof. +intros m e H. +unfold canonic. +now rewrite F2R_Zabs, canonic_exp_abs. +Qed. + +Theorem canonic_0: canonic (Float beta 0 (fexp (ln_beta beta 0%R))). +Proof. +unfold canonic; simpl; unfold canonic_exp. +replace (F2R {| Fnum := 0; Fexp := fexp (ln_beta beta 0) |}) with 0%R. +reflexivity. +unfold F2R; simpl; ring. +Qed. + + + Theorem canonic_unicity : forall f1 f2, canonic f1 -> @@ -756,6 +792,24 @@ refine (let H := _ in conj (proj1 H) (Rlt_le _ _ (proj2 H))). now apply mantissa_small_pos. Qed. + +Theorem exp_small_round_0_pos : + forall x ex, + (bpow (ex - 1) <= x < bpow ex)%R -> + round x =R0 -> (ex <= fexp ex)%Z . +Proof. +intros x ex H H1. +case (Zle_or_lt ex (fexp ex)); trivial; intros V. +contradict H1. +apply Rgt_not_eq. +apply Rlt_le_trans with (bpow (ex-1)). +apply bpow_gt_0. +apply (round_bounded_large_pos); assumption. +Qed. + + + + Theorem generic_format_round_pos : forall x, (0 < x)%R -> @@ -1014,6 +1068,24 @@ intros rnd' Hr x. apply round_bounded_large_pos... Qed. +Theorem exp_small_round_0 : + forall rnd {Hr : Valid_rnd rnd} x ex, + (bpow (ex - 1) <= Rabs x < bpow ex)%R -> + round rnd x =R0 -> (ex <= fexp ex)%Z . +Proof. +intros rnd Hr x ex H1 H2. +generalize Rabs_R0. +rewrite <- H2 at 1. +apply (round_abs_abs' (fun t rt => forall (ex : Z), +(bpow (ex - 1) <= t < bpow ex)%R -> +rt = 0%R -> (ex <= fexp ex)%Z)); trivial. +clear; intros rnd Hr x Hx. +now apply exp_small_round_0_pos. +Qed. + + + + Section monotone_abs. Variable rnd : R -> Z. @@ -1283,6 +1355,33 @@ rewrite <- mantissa_DN_small_pos with (1 := Hx) (2 := He). now rewrite <- canonic_exp_fexp_pos with (1 := Hx). Qed. + +Theorem round_DN_UP_lt : + forall x, ~ generic_format x -> + (round Zfloor x < x < round Zceil x)%R. +Proof with auto with typeclass_instances. +intros x Fx. +assert (Hx:(round Zfloor x <= x <= round Zceil x)%R). +split. +apply round_DN_pt. +apply round_UP_pt. +split. + destruct Hx as [Hx _]. + apply Rnot_le_lt; intro Hle. + assert (x = round Zfloor x) by now apply Rle_antisym. + rewrite H in Fx. + contradict Fx. + apply generic_format_round... +destruct Hx as [_ Hx]. +apply Rnot_le_lt; intro Hle. +assert (x = round Zceil x) by now apply Rle_antisym. +rewrite H in Fx. +contradict Fx. +apply generic_format_round... +Qed. + + + Theorem round_UP_small_pos : forall x ex, (bpow (ex - 1) <= x < bpow ex)%R -> diff --git a/flocq/Core/Fcore_rnd.v b/flocq/Core/Fcore_rnd.v index 6b4d807..94c9420 100644 --- a/flocq/Core/Fcore_rnd.v +++ b/flocq/Core/Fcore_rnd.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2010-2011 Sylvie Boldo +Copyright (C) 2010-2013 Sylvie Boldo #
# -Copyright (C) 2010-2011 Guillaume Melquiond +Copyright (C) 2010-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public diff --git a/flocq/Core/Fcore_rnd_ne.v b/flocq/Core/Fcore_rnd_ne.v index 0b0776e..6829c0c 100644 --- a/flocq/Core/Fcore_rnd_ne.v +++ b/flocq/Core/Fcore_rnd_ne.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2010-2011 Sylvie Boldo +Copyright (C) 2010-2013 Sylvie Boldo #
# -Copyright (C) 2010-2011 Guillaume Melquiond +Copyright (C) 2010-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public @@ -499,6 +499,25 @@ rewrite Zeven_plus. now rewrite eqb_sym. Qed. +Lemma round_NE_abs: + forall x : R, + round beta fexp ZnearestE (Rabs x) = Rabs (round beta fexp ZnearestE x). +Proof with auto with typeclass_instances. +intros x. +apply sym_eq. +unfold Rabs at 2. +destruct (Rcase_abs x) as [Hx|Hx]. +rewrite round_NE_opp. +apply Rabs_left1. +rewrite <- (round_0 beta fexp ZnearestE). +apply round_le... +now apply Rlt_le. +apply Rabs_pos_eq. +rewrite <- (round_0 beta fexp ZnearestE). +apply round_le... +now apply Rge_le. +Qed. + Theorem round_NE_pt : forall x, Rnd_NE_pt x (round beta fexp ZnearestE x). diff --git a/flocq/Core/Fcore_ulp.v b/flocq/Core/Fcore_ulp.v index 492fac6..07ef3ec 100644 --- a/flocq/Core/Fcore_ulp.v +++ b/flocq/Core/Fcore_ulp.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2010-2011 Sylvie Boldo +Copyright (C) 2010-2013 Sylvie Boldo #
# -Copyright (C) 2010-2011 Guillaume Melquiond +Copyright (C) 2010-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public @@ -1139,4 +1139,278 @@ apply le_pred_lt_aux ; try easy. now split. Qed. + +Theorem pred_succ : forall { monotone_exp : Monotone_exp fexp }, + forall x, F x -> (0 < x)%R -> pred (x + ulp x)=x. +Proof. +intros L x Fx Hx. +assert (x <= pred (x + ulp x))%R. +apply le_pred_lt. +assumption. +now apply generic_format_succ. +replace 0%R with (0+0)%R by ring; apply Rplus_lt_compat; try apply Hx. +apply bpow_gt_0. +apply Rplus_lt_reg_r with (-x)%R; ring_simplify. +apply bpow_gt_0. +apply Rle_antisym; trivial. +apply Rplus_le_reg_r with (ulp (pred (x + ulp x))). +rewrite pred_plus_ulp. +apply Rplus_le_compat_l. +now apply ulp_le. +replace 0%R with (0+0)%R by ring; apply Rplus_lt_compat; try apply Hx. +apply bpow_gt_0. +now apply generic_format_succ. +apply Rgt_not_eq. +now apply Rlt_le_trans with x. +Qed. + + +Theorem lt_UP_le_DN : + forall x y, F y -> + (y < round beta fexp Zceil x -> y <= round beta fexp Zfloor x)%R. +Proof with auto with typeclass_instances. +intros x y Fy Hlt. +apply round_DN_pt... +apply Rnot_lt_le. +contradict Hlt. +apply RIneq.Rle_not_lt. +apply round_UP_pt... +now apply Rlt_le. +Qed. + +Theorem pred_UP_le_DN : + forall x, (0 < round beta fexp Zceil x)%R -> + (pred (round beta fexp Zceil x) <= round beta fexp Zfloor x)%R. +Proof with auto with typeclass_instances. +intros x Pxu. +destruct (generic_format_EM beta fexp x) as [Fx|Fx]. +rewrite !round_generic... +now apply Rlt_le; apply pred_lt_id. +assert (let u := round beta fexp Zceil x in pred u < u)%R as Hup. + now apply pred_lt_id. +apply lt_UP_le_DN... +apply generic_format_pred... +now apply round_UP_pt. +Qed. + +Theorem pred_UP_eq_DN : + forall x, (0 < round beta fexp Zceil x)%R -> ~ F x -> + (pred (round beta fexp Zceil x) = round beta fexp Zfloor x)%R. +Proof with auto with typeclass_instances. +intros x Px Fx. +apply Rle_antisym. +now apply pred_UP_le_DN. +apply le_pred_lt; try apply generic_format_round... +pose proof round_DN_UP_lt _ _ _ Fx as HE. +now apply Rlt_trans with (1 := proj1 HE) (2 := proj2 HE). +Qed. + + + + + +(** Properties of rounding to nearest and ulp *) + +Theorem rnd_N_le_half_an_ulp: forall choice u v, + F u -> (0 < u)%R -> (v < u + (ulp u)/2)%R + -> (round beta fexp (Znearest choice) v <= u)%R. +Proof with auto with typeclass_instances. +intros choice u v Fu Hu H. +(* . *) +assert (0 < ulp u / 2)%R. +unfold Rdiv; apply Rmult_lt_0_compat. +unfold ulp; apply bpow_gt_0. +auto with real. +(* . *) +assert (ulp u / 2 < ulp u)%R. +apply Rlt_le_trans with (ulp u *1)%R;[idtac|right; ring]. +unfold Rdiv; apply Rmult_lt_compat_l. +apply bpow_gt_0. +apply Rmult_lt_reg_l with 2%R. +auto with real. +apply Rle_lt_trans with 1%R. +right; field. +rewrite Rmult_1_r; auto with real. +(* *) +apply Rnd_N_pt_monotone with F v (u + ulp u / 2)%R... +apply round_N_pt... +apply Rnd_DN_pt_N with (u+ulp u)%R. +pattern u at 3; replace u with (round beta fexp Zfloor (u + ulp u / 2)). +apply round_DN_pt... +apply round_DN_succ; try assumption. +split; try left; assumption. +replace (u+ulp u)%R with (round beta fexp Zceil (u + ulp u / 2)). +apply round_UP_pt... +apply round_UP_succ; try assumption... +split; try left; assumption. +right; field. +Qed. + + +Theorem rnd_N_ge_half_an_ulp_pred: forall choice u v, + F u -> (0 < pred u)%R -> (u - (ulp (pred u))/2 < v)%R + -> (u <= round beta fexp (Znearest choice) v)%R. +Proof with auto with typeclass_instances. +intros choice u v Fu Hu H. +(* . *) +assert (0 < u)%R. +apply Rlt_trans with (1:= Hu). +apply pred_lt_id. +assert (0 < ulp (pred u) / 2)%R. +unfold Rdiv; apply Rmult_lt_0_compat. +unfold ulp; apply bpow_gt_0. +auto with real. +assert (ulp (pred u) / 2 < ulp (pred u))%R. +apply Rlt_le_trans with (ulp (pred u) *1)%R;[idtac|right; ring]. +unfold Rdiv; apply Rmult_lt_compat_l. +apply bpow_gt_0. +apply Rmult_lt_reg_l with 2%R. +auto with real. +apply Rle_lt_trans with 1%R. +right; field. +rewrite Rmult_1_r; auto with real. +(* *) +apply Rnd_N_pt_monotone with F (u - ulp (pred u) / 2)%R v... +2: apply round_N_pt... +apply Rnd_UP_pt_N with (pred u). +pattern (pred u) at 2; replace (pred u) with (round beta fexp Zfloor (u - ulp (pred u) / 2)). +apply round_DN_pt... +replace (u - ulp (pred u) / 2)%R with (pred u + ulp (pred u) / 2)%R. +apply round_DN_succ; try assumption. +apply generic_format_pred; assumption. +split; [left|idtac]; assumption. +pattern u at 3; rewrite <- (pred_plus_ulp u); try assumption. +field. +now apply Rgt_not_eq. +pattern u at 3; replace u with (round beta fexp Zceil (u - ulp (pred u) / 2)). +apply round_UP_pt... +replace (u - ulp (pred u) / 2)%R with (pred u + ulp (pred u) / 2)%R. +apply trans_eq with (pred u +ulp(pred u))%R. +apply round_UP_succ; try assumption... +apply generic_format_pred; assumption. +split; [idtac|left]; assumption. +apply pred_plus_ulp; try assumption. +now apply Rgt_not_eq. +pattern u at 3; rewrite <- (pred_plus_ulp u); try assumption. +field. +now apply Rgt_not_eq. +pattern u at 4; rewrite <- (pred_plus_ulp u); try assumption. +right; field. +now apply Rgt_not_eq. +Qed. + + +Theorem rnd_N_ge_half_an_ulp: forall choice u v, + F u -> (0 < u)%R -> (u <> bpow (ln_beta beta u - 1))%R + -> (u - (ulp u)/2 < v)%R + -> (u <= round beta fexp (Znearest choice) v)%R. +Proof with auto with typeclass_instances. +intros choice u v Fu Hupos Hu H. +(* *) +assert (bpow (ln_beta beta u-1) <= pred u)%R. +apply le_pred_lt; try assumption. +apply generic_format_bpow. +assert (canonic_exp beta fexp u < ln_beta beta u)%Z. +apply ln_beta_generic_gt; try assumption. +now apply Rgt_not_eq. +unfold canonic_exp in H0. +ring_simplify (ln_beta beta u - 1 + 1)%Z. +omega. +destruct ln_beta as (e,He); simpl in *. +assert (bpow (e - 1) <= Rabs u)%R. +apply He. +now apply Rgt_not_eq. +rewrite Rabs_right in H0. +case H0; auto. +intros T; contradict T. +now apply sym_not_eq. +apply Rle_ge; now left. +assert (Hu2:(ulp (pred u) = ulp u)). +unfold ulp, canonic_exp. +apply f_equal; apply f_equal. +apply ln_beta_unique. +rewrite Rabs_right. +split. +assumption. +apply Rlt_trans with (1:=pred_lt_id _). +destruct ln_beta as (e,He); simpl in *. +rewrite Rabs_right in He. +apply He. +now apply Rgt_not_eq. +apply Rle_ge; now left. +apply Rle_ge, pred_ge_0; assumption. +apply rnd_N_ge_half_an_ulp_pred; try assumption. +apply Rlt_le_trans with (2:=H0). +apply bpow_gt_0. +rewrite Hu2; assumption. +Qed. + + +Lemma round_N_DN_betw: forall choice x d u, + Rnd_DN_pt (generic_format beta fexp) x d -> + Rnd_UP_pt (generic_format beta fexp) x u -> + (d<=x<(d+u)/2)%R -> + round beta fexp (Znearest choice) x = d. +Proof with auto with typeclass_instances. +intros choice x d u Hd Hu H. +apply Rnd_N_pt_unicity with (generic_format beta fexp) x d u; try assumption. +intros Y. +absurd (x < (d+u)/2)%R; try apply H. +apply Rle_not_lt; right. +apply Rplus_eq_reg_r with (-x)%R. +apply trans_eq with (- (x-d)/2 + (u-x)/2)%R. +field. +rewrite Y; field. +apply round_N_pt... +apply Rnd_DN_UP_pt_N with d u... +apply Hd. +right; apply trans_eq with (-(d-x))%R;[idtac|ring]. +apply Rabs_left1. +apply Rplus_le_reg_l with x; ring_simplify. +apply H. +rewrite Rabs_left1. +apply Rplus_le_reg_l with (d+x)%R. +apply Rmult_le_reg_l with (/2)%R. +auto with real. +apply Rle_trans with x. +right; field. +apply Rle_trans with ((d+u)/2)%R. +now left. +right; field. +apply Rplus_le_reg_l with x; ring_simplify. +apply H. +Qed. + + +Lemma round_N_UP_betw: forall choice x d u, + Rnd_DN_pt (generic_format beta fexp) x d -> + Rnd_UP_pt (generic_format beta fexp) x u -> + ((d+u)/2 < x <= u)%R -> + round beta fexp (Znearest choice) x = u. +Proof with auto with typeclass_instances. +intros choice x d u Hd Hu H. +rewrite <- (Ropp_involutive (round beta fexp (Znearest choice) x )), + <- (Ropp_involutive u) . +apply f_equal. +rewrite <- (Ropp_involutive x) . +rewrite round_N_opp, Ropp_involutive. +apply round_N_DN_betw with (-d)%R. +replace u with (round beta fexp Zceil x). +rewrite <- round_DN_opp. +apply round_DN_pt... +apply Rnd_UP_pt_unicity with (generic_format beta fexp) x... +apply round_UP_pt... +replace d with (round beta fexp Zfloor x). +rewrite <- round_UP_opp. +apply round_UP_pt... +apply Rnd_DN_pt_unicity with (generic_format beta fexp) x... +apply round_DN_pt... +split. +apply Ropp_le_contravar, H. +apply Rlt_le_trans with (-((d + u) / 2))%R. +apply Ropp_lt_contravar, H. +unfold Rdiv; right; ring. +Qed. + + End Fcore_ulp. diff --git a/flocq/Flocq_version.v b/flocq/Flocq_version.v index 662d83a..b375857 100644 --- a/flocq/Flocq_version.v +++ b/flocq/Flocq_version.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2011 Sylvie Boldo +Copyright (C) 2011-2013 Sylvie Boldo #
# -Copyright (C) 2011 Guillaume Melquiond +Copyright (C) 2011-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public @@ -28,4 +28,4 @@ Definition Flocq_version := Eval vm_compute in | String h t => parse t major (minor + N_of_ascii h - N_of_ascii "0"%char)%N | Empty_string => (major * 100 + minor)%N end in - parse "2.1.0"%string N0 N0. + parse "2.2.0"%string N0 N0. diff --git a/flocq/Prop/Fprop_Sterbenz.v b/flocq/Prop/Fprop_Sterbenz.v index 7d2c2e7..7260d2e 100644 --- a/flocq/Prop/Fprop_Sterbenz.v +++ b/flocq/Prop/Fprop_Sterbenz.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2010-2011 Sylvie Boldo +Copyright (C) 2010-2013 Sylvie Boldo #
# -Copyright (C) 2010-2011 Guillaume Melquiond +Copyright (C) 2010-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public diff --git a/flocq/Prop/Fprop_div_sqrt_error.v b/flocq/Prop/Fprop_div_sqrt_error.v index 84a0694..ec00ca4 100644 --- a/flocq/Prop/Fprop_div_sqrt_error.v +++ b/flocq/Prop/Fprop_div_sqrt_error.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2010-2011 Sylvie Boldo +Copyright (C) 2010-2013 Sylvie Boldo #
# -Copyright (C) 2010-2011 Guillaume Melquiond +Copyright (C) 2010-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public diff --git a/flocq/Prop/Fprop_mult_error.v b/flocq/Prop/Fprop_mult_error.v index e3e094c..e84e80b 100644 --- a/flocq/Prop/Fprop_mult_error.v +++ b/flocq/Prop/Fprop_mult_error.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2010-2011 Sylvie Boldo +Copyright (C) 2010-2013 Sylvie Boldo #
# -Copyright (C) 2010-2011 Guillaume Melquiond +Copyright (C) 2010-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public diff --git a/flocq/Prop/Fprop_plus_error.v b/flocq/Prop/Fprop_plus_error.v index d9dee7c..ddae698 100644 --- a/flocq/Prop/Fprop_plus_error.v +++ b/flocq/Prop/Fprop_plus_error.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2010-2011 Sylvie Boldo +Copyright (C) 2010-2013 Sylvie Boldo #
# -Copyright (C) 2010-2011 Guillaume Melquiond +Copyright (C) 2010-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public diff --git a/flocq/Prop/Fprop_relative.v b/flocq/Prop/Fprop_relative.v index 8df7336..a8cc1ff 100644 --- a/flocq/Prop/Fprop_relative.v +++ b/flocq/Prop/Fprop_relative.v @@ -2,9 +2,9 @@ This file is part of the Flocq formalization of floating-point arithmetic in Coq: http://flocq.gforge.inria.fr/ -Copyright (C) 2010-2011 Sylvie Boldo +Copyright (C) 2010-2013 Sylvie Boldo #
# -Copyright (C) 2010-2011 Guillaume Melquiond +Copyright (C) 2010-2013 Guillaume Melquiond This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public @@ -35,9 +35,9 @@ Section relative_error_conversion. Variable rnd : R -> Z. Context { valid_rnd : Valid_rnd rnd }. -Lemma relative_error_lt_conversion : +Lemma relative_error_lt_conversion' : forall x b, (0 < b)%R -> - (Rabs (round beta fexp rnd x - x) < b * Rabs x)%R -> + (x <> 0 -> Rabs (round beta fexp rnd x - x) < b * Rabs x)%R -> exists eps, (Rabs eps < b)%R /\ round beta fexp rnd x = (x * (1 + eps))%R. Proof with auto with typeclass_instances. @@ -50,6 +50,7 @@ now rewrite Rabs_R0. rewrite Hx0, Rmult_0_l. apply round_0... (* *) +specialize (Hxb Hx0). exists ((round beta fexp rnd x - x) / x)%R. split. 2: now field. unfold Rdiv. @@ -61,6 +62,19 @@ rewrite Rinv_l with (1 := Hx0). now rewrite Rabs_R1, Rmult_1_r. Qed. +(* TODO: remove *) +Lemma relative_error_lt_conversion : + forall x b, (0 < b)%R -> + (Rabs (round beta fexp rnd x - x) < b * Rabs x)%R -> + exists eps, + (Rabs eps < b)%R /\ round beta fexp rnd x = (x * (1 + eps))%R. +Proof. +intros x b Hb0 Hxb. +apply relative_error_lt_conversion'. +exact Hb0. +now intros _. +Qed. + Lemma relative_error_le_conversion : forall x b, (0 <= b)%R -> (Rabs (round beta fexp rnd x - x) <= b * Rabs x)%R -> @@ -154,18 +168,28 @@ rewrite F2R_0, F2R_Zabs. now apply Rabs_pos_lt. Qed. -Theorem relative_error_F2R_emin_ex : +Theorem relative_error_F2R_emin_ex' : forall m, let x := F2R (Float beta m emin) in - (x <> 0)%R -> exists eps, (Rabs eps < bpow (-p + 1))%R /\ round beta fexp rnd x = (x * (1 + eps))%R. Proof with auto with typeclass_instances. -intros m x Hx. -apply relative_error_lt_conversion... +intros m x. +apply relative_error_lt_conversion'... apply bpow_gt_0. now apply relative_error_F2R_emin. Qed. +(* TODO: remove *) +Theorem relative_error_F2R_emin_ex : + forall m, let x := F2R (Float beta m emin) in + (x <> 0)%R -> + exists eps, + (Rabs eps < bpow (-p + 1))%R /\ round beta fexp rnd x = (x * (1 + eps))%R. +Proof with auto with typeclass_instances. +intros m x _. +apply relative_error_F2R_emin_ex'. +Qed. + Theorem relative_error_round : (0 < p)%Z -> forall x, @@ -404,6 +428,7 @@ Qed. Variable rnd : R -> Z. Context { valid_rnd : Valid_rnd rnd }. +(* TODO: remove *) Theorem relative_error_FLT_F2R_emin : forall m, let x := F2R (Float beta m (emin + prec - 1)) in (x <> 0)%R -> @@ -424,14 +449,49 @@ apply relative_error with (emin + prec - 1)%Z... apply relative_error_FLT_aux. Qed. +Theorem relative_error_FLT_F2R_emin' : + forall m, let x := F2R (Float beta m emin) in + (x <> 0)%R -> + (Rabs (round beta (FLT_exp emin prec) rnd x - x) < bpow (-prec + 1) * Rabs x)%R. +Proof with auto with typeclass_instances. +intros m x Zx. +destruct (Rlt_or_le (Rabs x) (bpow (emin + prec - 1))) as [Hx|Hx]. +rewrite round_generic... +unfold Rminus. +rewrite Rplus_opp_r, Rabs_R0. +apply Rmult_lt_0_compat. +apply bpow_gt_0. +now apply Rabs_pos_lt. +apply generic_format_FLT_FIX... +apply Rlt_le. +apply Rlt_le_trans with (1 := Hx). +apply bpow_le. +apply Zle_pred. +apply generic_format_FIX. +now exists (Float beta m emin). +now apply relative_error_FLT. +Qed. + +Theorem relative_error_FLT_F2R_emin_ex' : + forall m, let x := F2R (Float beta m emin) in + exists eps, + (Rabs eps < bpow (-prec + 1))%R /\ round beta (FLT_exp emin prec) rnd x = (x * (1 + eps))%R. +Proof with auto with typeclass_instances. +intros m x. +apply relative_error_lt_conversion'... +apply bpow_gt_0. +now apply relative_error_FLT_F2R_emin'. +Qed. + +(* TODO: remove *) Theorem relative_error_FLT_F2R_emin_ex : forall m, let x := F2R (Float beta m (emin + prec - 1)) in (x <> 0)%R -> exists eps, (Rabs eps < bpow (-prec + 1))%R /\ round beta (FLT_exp emin prec) rnd x = (x * (1 + eps))%R. Proof with auto with typeclass_instances. -intros m x Hx. -apply relative_error_lt_conversion... +intros m x _. +apply relative_error_lt_conversion'... apply bpow_gt_0. now apply relative_error_FLT_F2R_emin. Qed. @@ -488,6 +548,32 @@ apply relative_error_N_round with (emin + prec - 1)%Z... apply relative_error_FLT_aux. Qed. +Theorem relative_error_N_FLT_F2R_emin' : + forall m, let x := F2R (Float beta m emin) in + (Rabs (round beta (FLT_exp emin prec) (Znearest choice) x - x) <= /2 * bpow (-prec + 1) * Rabs x)%R. +Proof with auto with typeclass_instances. +intros m x. +destruct (Rlt_or_le (Rabs x) (bpow (emin + prec - 1))) as [Hx|Hx]. +rewrite round_generic... +unfold Rminus. +rewrite Rplus_opp_r, Rabs_R0. +apply Rmult_le_pos. +apply Rmult_le_pos. +apply Rlt_le. +apply (RinvN_pos 1). +apply bpow_ge_0. +apply Rabs_pos. +apply generic_format_FLT_FIX... +apply Rlt_le. +apply Rlt_le_trans with (1 := Hx). +apply bpow_le. +apply Zle_pred. +apply generic_format_FIX. +now exists (Float beta m emin). +now apply relative_error_N_FLT. +Qed. + +(* TODO: remove *) Theorem relative_error_N_FLT_F2R_emin : forall m, let x := F2R (Float beta m (emin + prec - 1)) in (Rabs (round beta (FLT_exp emin prec) (Znearest choice) x - x) <= /2 * bpow (-prec + 1) * Rabs x)%R. @@ -497,6 +583,21 @@ apply relative_error_N_F2R_emin... apply relative_error_FLT_aux. Qed. +Theorem relative_error_N_FLT_F2R_emin_ex' : + forall m, let x := F2R (Float beta m emin) in + exists eps, + (Rabs eps <= /2 * bpow (-prec + 1))%R /\ round beta (FLT_exp emin prec) (Znearest choice) x = (x * (1 + eps))%R. +Proof with auto with typeclass_instances. +intros m x. +apply relative_error_le_conversion... +apply Rmult_le_pos. +apply Rlt_le. +apply (RinvN_pos 1). +apply bpow_ge_0. +now apply relative_error_N_FLT_F2R_emin'. +Qed. + +(* TODO: remove *) Theorem relative_error_N_FLT_F2R_emin_ex : forall m, let x := F2R (Float beta m (emin + prec - 1)) in exists eps, @@ -512,6 +613,33 @@ apply bpow_gt_0. now apply relative_error_N_FLT_F2R_emin. Qed. +Theorem relative_error_N_FLT_round_F2R_emin' : + forall m, let x := F2R (Float beta m emin) in + (Rabs (round beta (FLT_exp emin prec) (Znearest choice) x - x) <= /2 * bpow (-prec + 1) * Rabs (round beta (FLT_exp emin prec) (Znearest choice) x))%R. +Proof with auto with typeclass_instances. +intros m x. +destruct (Rlt_or_le (Rabs x) (bpow (emin + prec - 1))) as [Hx|Hx]. +rewrite round_generic... +unfold Rminus. +rewrite Rplus_opp_r, Rabs_R0. +apply Rmult_le_pos. +apply Rmult_le_pos. +apply Rlt_le. +apply (RinvN_pos 1). +apply bpow_ge_0. +apply Rabs_pos. +apply generic_format_FLT_FIX... +apply Rlt_le. +apply Rlt_le_trans with (1 := Hx). +apply bpow_le. +apply Zle_pred. +apply generic_format_FIX. +now exists (Float beta m emin). +apply relative_error_N_round with (emin := (emin + prec - 1)%Z)... +apply relative_error_FLT_aux. +Qed. + +(* TODO: remove *) Theorem relative_error_N_FLT_round_F2R_emin : forall m, let x := F2R (Float beta m (emin + prec - 1)) in (Rabs (round beta (FLT_exp emin prec) (Znearest choice) x - x) <= /2 * bpow (-prec + 1) * Rabs (round beta (FLT_exp emin prec) (Znearest choice) x))%R. @@ -606,18 +734,28 @@ apply He. Qed. (** 1+#ε# property in any rounding in FLX *) -Theorem relative_error_FLX_ex : +Theorem relative_error_FLX_ex' : forall x, - (x <> 0)%R -> exists eps, (Rabs eps < bpow (-prec + 1))%R /\ round beta (FLX_exp prec) rnd x = (x * (1 + eps))%R. Proof with auto with typeclass_instances. -intros x Hx. -apply relative_error_lt_conversion... +intros x. +apply relative_error_lt_conversion'... apply bpow_gt_0. now apply relative_error_FLX. Qed. +(* TODO: remove *) +Theorem relative_error_FLX_ex : + forall x, + (x <> 0)%R -> + exists eps, + (Rabs eps < bpow (-prec + 1))%R /\ round beta (FLX_exp prec) rnd x = (x * (1 + eps))%R. +Proof with auto with typeclass_instances. +intros x _. +apply relative_error_FLX_ex'. +Qed. + Theorem relative_error_FLX_round : forall x, (x <> 0)%R -> diff --git a/ia32/Nan.v b/ia32/Nan.v new file mode 100644 index 0000000..dadf0ca --- /dev/null +++ b/ia32/Nan.v @@ -0,0 +1,28 @@ +Require Import Fappli_IEEE. +Require Import Fappli_IEEE_bits. +Require Import Floats. +Require Import ZArith. +Require Import Integers. + +(* Needed to break a circular reference after extraction *) +Definition transform_quiet_pl := + Eval unfold Float.transform_quiet_pl in Float.transform_quiet_pl. + +Program Definition default_pl : bool * nan_pl 53 := (true, nat_iter 51 xO xH). + +Definition binop_pl (pl1 pl2:binary64) : bool*nan_pl 53 := + match pl1, pl2 with + | B754_nan s1 pl1, _ => (s1, transform_quiet_pl pl1) + | _, B754_nan s2 pl2 => (s2, transform_quiet_pl pl2) + | _, _ => default_pl + end. + +Theorem binop_propagate1: Float.binop_propagate1_prop binop_pl. +Proof. + repeat intro. destruct f1, f2; try discriminate; simpl; reflexivity. +Qed. + +Theorem binop_propagate2: Float.binop_propagate2_prop binop_pl. +Proof. + repeat intro. destruct f1, f2, f3; try discriminate; simpl; reflexivity. +Qed. diff --git a/lib/Floats.v b/lib/Floats.v index 94c19d2..0151cf0 100644 --- a/lib/Floats.v +++ b/lib/Floats.v @@ -38,20 +38,69 @@ Definition zero: float := B754_zero _ _ false. (**r the float [+0.0] *) Definition eq_dec: forall (f1 f2: float), {f1 = f2} + {f1 <> f2}. Proof. Ltac try_not_eq := try solve [right; congruence]. - destruct f1, f2; + destruct f1 as [| |? []|], f2 as [| |? []|]; try destruct b; try destruct b0; - try solve [left; auto]; try_not_eq; + try solve [left; auto]; try_not_eq. + destruct (positive_eq_dec x x0); try_not_eq; + subst; left; f_equal; f_equal; apply proof_irr. + destruct (positive_eq_dec x x0); try_not_eq; + subst; left; f_equal; f_equal; apply proof_irr. + destruct (positive_eq_dec m m0); try_not_eq; + destruct (Z_eq_dec e e1); try solve [right; intro H; inv H; congruence]; + subst; left; rewrite (proof_irr e0 e2); auto. destruct (positive_eq_dec m m0); try_not_eq; destruct (Z_eq_dec e e1); try solve [right; intro H; inv H; congruence]; subst; left; rewrite (proof_irr e0 e2); auto. Defined. +(* Transform a Nan payload to a quiet Nan payload. + This is not part of the IEEE754 standard, but shared between all + architectures of Compcert. *) +Program Definition transform_quiet_pl (pl:nan_pl 53) : nan_pl 53 := + Pos.lor pl (nat_iter 51 xO xH). +Next Obligation. + destruct pl. + simpl. rewrite Z.ltb_lt in *. + assert (forall x, S (Fcore_digits.digits2_Pnat x) = Pos.to_nat (Pos.size x)). + { induction x0; simpl; auto; rewrite IHx0; zify; omega. } + fold (Z.of_nat (S (Fcore_digits.digits2_Pnat (Pos.lor x 2251799813685248)))). + rewrite H, positive_nat_Z, Psize_log_inf, <- Zlog2_log_inf in *. clear H. + change (Z.pos (Pos.lor x 2251799813685248)) with (Z.lor (Z.pos x) 2251799813685248%Z). + rewrite Z.log2_lor by (zify; omega). + apply Z.max_case. auto. simpl. omega. +Qed. + +Lemma nan_payload_fequal: + forall prec p1 e1 p2 e2, p1 = p2 -> (exist _ p1 e1:nan_pl prec) = exist _ p2 e2. +Proof. + simpl; intros; subst. f_equal. apply Fcore_Zaux.eqbool_irrelevance. +Qed. + +Lemma lor_idempotent: + forall x y, Pos.lor (Pos.lor x y) y = Pos.lor x y. +Proof. + induction x; destruct y; simpl; f_equal; auto; + induction y; simpl; f_equal; auto. +Qed. + +Lemma transform_quiet_pl_idempotent: + forall pl, transform_quiet_pl (transform_quiet_pl pl) = transform_quiet_pl pl. +Proof. + intros []; simpl; intros. apply nan_payload_fequal. + simpl. apply lor_idempotent. +Qed. + (** Arithmetic operations *) -Definition neg: float -> float := b64_opp. (**r opposite (change sign) *) +(* The Nan payload operations for neg and abs is not part of the IEEE754 + standard, but shared between all architectures of Compcert. *) +Definition neg_pl (s:bool) (pl:nan_pl 53) := (negb s, pl). +Definition abs_pl (s:bool) (pl:nan_pl 53) := (false, pl). + +Definition neg: float -> float := b64_opp neg_pl. (**r opposite (change sign) *) Definition abs (x: float): float := (**r absolute value (set sign to [+]) *) match x with - | B754_nan => x + | B754_nan s pl => let '(s, pl) := abs_pl s pl in B754_nan _ _ s pl | B754_infinity _ => B754_infinity _ _ false | B754_finite _ m e H => B754_finite _ _ false m e H | B754_zero _ => B754_zero _ _ false @@ -71,9 +120,30 @@ Definition binary_normalize32_correct (m e:Z) (s:bool) := binary_normalize_correct 24 128 eq_refl eq_refl mode_NE m e s. Global Opaque binary_normalize32_correct. +(* The Nan payload operations for single <-> double conversions are not part of + the IEEE754 standard, but shared between all architectures of Compcert. *) +Definition floatofbinary32_pl (s:bool) (pl:nan_pl 24) : (bool * nan_pl 53). + refine (s, transform_quiet_pl (exist _ (Pos.shiftl_nat (proj1_sig pl) 29) _)). + abstract ( + destruct pl; unfold proj1_sig, Pos.shiftl_nat, nat_iter, Fcore_digits.digits2_Pnat; + fold (Fcore_digits.digits2_Pnat x); + rewrite Z.ltb_lt in *; + zify; omega). +Defined. + +Definition binary32offloat_pl (s:bool) (pl:nan_pl 53) : (bool * nan_pl 24). + refine (s, exist _ (Pos.shiftr_nat (proj1_sig (transform_quiet_pl pl)) 29) _). + abstract ( + destruct (transform_quiet_pl pl); unfold proj1_sig, Pos.shiftr_nat, nat_iter; + rewrite Z.ltb_lt in *; + assert (forall x, Fcore_digits.digits2_Pnat (Pos.div2 x) = + (Fcore_digits.digits2_Pnat x - 1)%nat) by (destruct x0; simpl; zify; omega); + rewrite !H, <- !NPeano.Nat.sub_add_distr; zify; omega). +Defined. + Definition floatofbinary32 (f: binary32) : float := (**r single precision embedding in double precision *) match f with - | B754_nan => B754_nan _ _ + | B754_nan s pl => let '(s, pl) := floatofbinary32_pl s pl in B754_nan _ _ s pl | B754_infinity s => B754_infinity _ _ s | B754_zero s => B754_zero _ _ s | B754_finite s m e _ => @@ -82,7 +152,7 @@ Definition floatofbinary32 (f: binary32) : float := (**r single precision embedd Definition binary32offloat (f: float) : binary32 := (**r conversion to single precision *) match f with - | B754_nan => B754_nan _ _ + | B754_nan s pl => let '(s, pl) := binary32offloat_pl s pl in B754_nan _ _ s pl | B754_infinity s => B754_infinity _ _ s | B754_zero s => B754_zero _ _ s | B754_finite s m e _ => @@ -190,15 +260,33 @@ Definition singleoflong (n:int64): float := (**r conversion from signed 64-bit i Definition singleoflongu (n:int64): float:= (**r conversion from unsigned 64-bit int to single-precision float *) floatofbinary32 (binary_normalize32 (Int64.unsigned n) 0 false). -Definition add: float -> float -> float := b64_plus mode_NE. (**r addition *) -Definition sub: float -> float -> float := b64_minus mode_NE. (**r subtraction *) -Definition mul: float -> float -> float := b64_mult mode_NE. (**r multiplication *) -Definition div: float -> float -> float := b64_div mode_NE. (**r division *) +Parameter binop_pl: binary64 -> binary64 -> bool*nan_pl 53. + +Definition binop_propagate1_prop binop_pl := + forall f1 f2:float, is_nan _ _ f1 = true -> is_nan _ _ f2 = false -> + binop_pl f1 f1 = (binop_pl f1 f2 : bool*nan_pl 53). + + +(* Proved in Nan.v for different architectures *) +Hypothesis binop_propagate1: binop_propagate1_prop binop_pl. + +Definition binop_propagate2_prop binop_pl := + forall f1 f2 f3:float, is_nan _ _ f1 = true -> + is_nan _ _ f2 = false -> is_nan _ _ f3 = false -> + binop_pl f1 f2 = (binop_pl f1 f3 : bool*nan_pl 53). + +(* Proved in Nan.v for different architectures *) +Hypothesis binop_propagate2: binop_propagate2_prop binop_pl. + +Definition add: float -> float -> float := b64_plus binop_pl mode_NE. (**r addition *) +Definition sub: float -> float -> float := b64_minus binop_pl mode_NE. (**r subtraction *) +Definition mul: float -> float -> float := b64_mult binop_pl mode_NE. (**r multiplication *) +Definition div: float -> float -> float := b64_div binop_pl mode_NE. (**r division *) Definition order_float (f1 f2:float): option Datatypes.comparison := match f1, f2 with - | B754_nan,_ | _,B754_nan => None - | B754_infinity true, B754_infinity true + | B754_nan _ _,_ | _,B754_nan _ _ => None + | B754_infinity true, B754_infinity true | B754_infinity false, B754_infinity false => Some Eq | B754_infinity true, _ => Some Lt | B754_infinity false, _ => Some Gt @@ -229,7 +317,7 @@ Definition order_float (f1 f2:float): option Datatypes.comparison := end. Definition cmp (c:comparison) (f1 f2:float) : bool := (**r comparison *) - match c with + match c with | Ceq => match order_float f1 f2 with Some Eq => true | _ => false end | Cne => @@ -242,7 +330,7 @@ Definition cmp (c:comparison) (f1 f2:float) : bool := (**r comparison *) match order_float f1 f2 with Some Gt => true | _ => false end | Cge => match order_float f1 f2 with Some(Gt|Eq) => true | _ => false end - end. + end. (** Conversions between floats and their concrete in-memory representation as a sequence of 64 bits (double precision) or 32 bits (single precision). *) @@ -264,29 +352,17 @@ Ltac compute_this val := let x := fresh in set val as x in *; vm_compute in x; subst x. Ltac smart_omega := - simpl radix_val in *; simpl Zpower in *; + simpl radix_val in *; simpl Zpower in *; compute_this Int.modulus; compute_this Int.half_modulus; compute_this Int.max_unsigned; compute_this Int.min_signed; compute_this Int.max_signed; compute_this Int64.modulus; compute_this Int64.half_modulus; compute_this Int64.max_unsigned; compute_this (Zpower_pos 2 1024); compute_this (Zpower_pos 2 53); compute_this (Zpower_pos 2 52); - omega. - -Theorem addf_commut: forall f1 f2, add f1 f2 = add f2 f1. -Proof. - intros. - destruct f1, f2; simpl; try reflexivity; try (destruct b, b0; reflexivity). - rewrite Zplus_comm; rewrite Zmin_comm; reflexivity. -Qed. - -Theorem subf_addf_opp: forall f1 f2, sub f1 f2 = add f1 (neg f2). -Proof. - destruct f1, f2; reflexivity. -Qed. + zify; omega. Lemma floatofbinary32_exact : - forall f, is_finite_strict _ _ f = true -> + forall f, is_finite_strict _ _ f = true -> is_finite_strict _ _ (floatofbinary32 f) = true /\ B2R _ _ f = B2R _ _ (floatofbinary32 f). Proof. destruct f as [ | | |s m e]; try discriminate; intro. @@ -305,10 +381,12 @@ Proof. now apply bpow_lt. Qed. -Lemma binary32offloatofbinary32 : - forall f, binary32offloat (floatofbinary32 f) = f. +Lemma binary32offloatofbinary32_num : + forall f, is_nan _ _ f = false -> + binary32offloat (floatofbinary32 f) = f. Proof. - intro; pose proof (floatofbinary32_exact f); destruct f as [ | | |s m e]; try reflexivity. + intros f Hnan; pose proof (floatofbinary32_exact f); destruct f as [ | | |s m e]; try reflexivity. + discriminate. specialize (H eq_refl); destruct H. destruct (floatofbinary32 (B754_finite 24 128 s m e e0)) as [ | | |s1 m1 e1]; try discriminate. unfold binary32offloat. @@ -324,22 +402,76 @@ Proof. now apply generic_format_B2R. Qed. +Lemma floatofbinary32offloatofbinary32_pl: + forall s pl, + prod_rect (fun _ => _) floatofbinary32_pl (prod_rect (fun _ => _) binary32offloat_pl (floatofbinary32_pl s pl)) = floatofbinary32_pl s pl. +Proof. + destruct pl. unfold binary32offloat_pl, floatofbinary32_pl. + unfold transform_quiet_pl, proj1_sig. simpl. + f_equal. apply nan_payload_fequal. + unfold Pos.shiftr_nat. simpl. + rewrite !lor_idempotent. reflexivity. +Qed. + +Lemma floatofbinary32offloatofbinary32 : + forall f, floatofbinary32 (binary32offloat (floatofbinary32 f)) = floatofbinary32 f. +Proof. + destruct f; try (rewrite binary32offloatofbinary32_num; tauto). + unfold floatofbinary32, binary32offloat. + rewrite <- floatofbinary32offloatofbinary32_pl at 2. + reflexivity. +Qed. + +Lemma binary32offloatofbinary32offloat_pl: + forall s pl, + prod_rect (fun _ => _) binary32offloat_pl (prod_rect (fun _ => _) floatofbinary32_pl (binary32offloat_pl s pl)) = binary32offloat_pl s pl. +Proof. + destruct pl. unfold binary32offloat_pl, floatofbinary32_pl. unfold prod_rect. + f_equal. apply nan_payload_fequal. + rewrite transform_quiet_pl_idempotent. + unfold transform_quiet_pl, proj1_sig. + change 51 with (29+22). + clear - x. revert x. unfold Pos.shiftr_nat, Pos.shiftl_nat. + induction (29)%nat. intro. simpl. apply lor_idempotent. + intro. + rewrite !nat_iter_succ_r with (f:=Pos.div2). + destruct x; simpl; try apply IHn. + clear IHn. induction n. reflexivity. + rewrite !nat_iter_succ_r with (f:=Pos.div2). auto. +Qed. + +Lemma binary32offloatofbinary32offloat : + forall f, binary32offloat (floatofbinary32 (binary32offloat f)) = binary32offloat f. +Proof. + destruct f; try (rewrite binary32offloatofbinary32_num; simpl; tauto). + unfold floatofbinary32, binary32offloat. + rewrite <- binary32offloatofbinary32offloat_pl at 2. + reflexivity. + rewrite binary32offloatofbinary32_num; simpl. auto. + unfold binary_normalize32. + pose proof (binary_normalize32_correct (cond_Zopp b (Z.pos m)) e b). + destruct binary_normalize; auto. simpl in H. + destruct Rlt_bool in H. intuition. + unfold binary_overflow in H. destruct n. + destruct overflow_to_inf in H; discriminate. +Qed. + Theorem singleoffloat_idem: forall f, singleoffloat (singleoffloat f) = singleoffloat f. Proof. - intros; unfold singleoffloat; rewrite binary32offloatofbinary32; reflexivity. + intros; unfold singleoffloat; rewrite binary32offloatofbinary32offloat; reflexivity. Qed. Theorem singleoflong_idem: forall n, singleoffloat (singleoflong n) = singleoflong n. Proof. - intros; unfold singleoffloat, singleoflong; rewrite binary32offloatofbinary32; reflexivity. + intros; unfold singleoffloat, singleoflong. rewrite floatofbinary32offloatofbinary32; reflexivity. Qed. Theorem singleoflongu_idem: forall n, singleoffloat (singleoflongu n) = singleoflongu n. Proof. - intros; unfold singleoffloat, singleoflongu; rewrite binary32offloatofbinary32; reflexivity. + intros; unfold singleoffloat, singleoflongu. rewrite floatofbinary32offloatofbinary32; reflexivity. Qed. Definition is_single (f: float) : Prop := exists s, f = floatofbinary32 s. @@ -353,8 +485,8 @@ Qed. Theorem singleoffloat_of_single: forall f, is_single f -> singleoffloat f = f. Proof. - intros. destruct H as [s EQ]. subst f. unfold singleoffloat. - rewrite binary32offloatofbinary32; reflexivity. + intros. destruct H as [s EQ]. subst f. unfold singleoffloat. + apply floatofbinary32offloatofbinary32. Qed. Theorem is_single_dec: forall f, {is_single f} + {~is_single f}. @@ -402,9 +534,9 @@ Proof. replace 0%R with (0*bpow radix2 e0)%R by ring; apply Rmult_lt_compat_r; [apply bpow_gt_0; reflexivity|now apply (Z2R_lt _ 0)]. apply Rmult_gt_0_compat; [apply (Z2R_lt 0); reflexivity|now apply bpow_gt_0]. - destruct b, b0; try (now apply_Rcompare; apply H5); inversion H3; + destruct b, b0; try (now apply_Rcompare; apply H5); inversion H3; try (apply_Rcompare; apply H4; rewrite H, H1 in H7; assumption); - try (apply_Rcompare; do 2 rewrite Z2R_opp, Ropp_mult_distr_l_reverse; + try (apply_Rcompare; do 2 rewrite Z2R_opp, Ropp_mult_distr_l_reverse; apply Ropp_lt_contravar; apply H4; rewrite H, H1 in H7; assumption); rewrite H7, Rcompare_mult_r, Rcompare_Z2R by (apply bpow_gt_0); reflexivity. Qed. @@ -472,7 +604,16 @@ Proof. destruct f. simpl; try destruct b; vm_compute; split; congruence. simpl; try destruct b; vm_compute; split; congruence. - simpl; vm_compute; split; congruence. + destruct n as [p Hp]. + simpl. rewrite Z.ltb_lt in Hp. + apply Zlt_succ_le with (m:=52) in Hp. + apply Zpower_le with (r:=radix2) in Hp. + edestruct Fcore_digits.digits2_Pnat_correct. + rewrite Zpower_nat_Z in H0. + eapply Z.lt_le_trans in Hp; eauto. + unfold join_bits; destruct b. + compute_this ((2 ^ 11 + 2047) * 2 ^ 52). smart_omega. + compute_this ((0 + 2047) * 2 ^ 52). smart_omega. unfold bits_of_binary_float, join_bits. destruct (andb_prop _ _ e0); apply Zle_bool_imp_le in H0; apply Zeq_bool_eq in H; unfold FLT_exp in H. match goal with [H:Zmax ?x ?y = e|-_] => pose proof (Zle_max_l x y); pose proof (Zle_max_r x y) end. @@ -494,7 +635,17 @@ Proof. destruct (binary32offloat f). simpl; try destruct b; vm_compute; split; congruence. simpl; try destruct b; vm_compute; split; congruence. - simpl; vm_compute; split; congruence. + destruct n as [p Hp]. + simpl. rewrite Z.ltb_lt in Hp. + apply Zlt_succ_le with (m:=23) in Hp. + apply Zpower_le with (r:=radix2) in Hp. + edestruct Fcore_digits.digits2_Pnat_correct. + rewrite Zpower_nat_Z in H0. + eapply Z.lt_le_trans in Hp; eauto. + compute_this (radix2^23). + unfold join_bits; destruct b. + compute_this ((2 ^ 8 + 255) * 2 ^ 23). smart_omega. + compute_this ((0 + 255) * 2 ^ 23). smart_omega. unfold bits_of_binary_float, join_bits. destruct (andb_prop _ _ e0); apply Zle_bool_imp_le in H0; apply Zeq_bool_eq in H. unfold FLT_exp in H. @@ -512,13 +663,13 @@ Qed. Theorem bits_of_singleoffloat: forall f, bits_of_single (singleoffloat f) = bits_of_single f. Proof. - intro; unfold singleoffloat, bits_of_single; rewrite binary32offloatofbinary32; reflexivity. + intro; unfold singleoffloat, bits_of_single; rewrite binary32offloatofbinary32offloat; reflexivity. Qed. Theorem singleoffloat_of_bits: forall b, singleoffloat (single_of_bits b) = single_of_bits b. Proof. - intro; unfold singleoffloat, single_of_bits; rewrite binary32offloatofbinary32; reflexivity. + intro; unfold singleoffloat, single_of_bits; rewrite floatofbinary32offloatofbinary32; reflexivity. Qed. Theorem single_of_bits_is_single: @@ -539,7 +690,7 @@ Lemma round_exact: round radix2 (FLT_exp (3 - 1024 - 53) 53) (round_mode mode_NE) (Z2R n) = Z2R n. Proof. - intros; rewrite round_generic; [reflexivity|now apply valid_rnd_round_mode|]. + intros; rewrite round_generic; [reflexivity|now apply valid_rnd_round_mode|]. apply generic_format_FLT; exists (Float radix2 n 0). unfold F2R, Fnum, Fexp, bpow; rewrite Rmult_1_r; intuition. pose proof (Zabs_spec n); now smart_omega. @@ -551,7 +702,7 @@ Lemma binary_normalize64_exact: is_finite _ _ (binary_normalize64 n 0 false) = true. Proof. intros; pose proof (binary_normalize64_correct n 0 false). - unfold F2R, Fnum, Fexp, bpow in H0; rewrite Rmult_1_r, round_exact, Rlt_bool_true in H0; try assumption. + unfold F2R, Fnum, Fexp, bpow in H0; rewrite Rmult_1_r, round_exact, Rlt_bool_true in H0; try now intuition. rewrite <- Z2R_abs; apply Z2R_lt; pose proof (Zabs_spec n); now smart_omega. Qed. @@ -560,7 +711,7 @@ Theorem floatofintu_floatofint_1: Int.ltu x ox8000_0000 = true -> floatofintu x = floatofint x. Proof. - unfold floatofintu, floatofint, Int.signed, Int.ltu; intro. + unfold floatofintu, floatofint, Int.signed, Int.ltu; intro. change (Int.unsigned ox8000_0000) with Int.half_modulus. destruct (zlt (Int.unsigned x) Int.half_modulus); now intuition. Qed. @@ -571,19 +722,19 @@ Theorem floatofintu_floatofint_2: floatofintu x = add (floatofint (Int.sub x ox8000_0000)) (floatofintu ox8000_0000). Proof. - unfold floatofintu, floatofint, Int.signed, Int.ltu, Int.sub; intros. + unfold floatofintu, floatofint, Int.signed, Int.ltu, Int.sub; intros. pose proof (Int.unsigned_range x). compute_this (Int.unsigned ox8000_0000). destruct (zlt (Int.unsigned x) 2147483648); try discriminate. rewrite Int.unsigned_repr by smart_omega. destruct (zlt ((Int.unsigned x) - 2147483648) Int.half_modulus). - unfold add, b64_plus. - match goal with [|- _ = Bplus _ _ _ _ _ ?x ?y] => - pose proof (Bplus_correct 53 1024 eq_refl eq_refl mode_NE x y) end. + unfold add, b64_plus. + match goal with [|- _ = Bplus _ _ _ _ _ _ ?x ?y] => + pose proof (Bplus_correct 53 1024 eq_refl eq_refl binop_pl mode_NE x y) end. do 2 rewrite (fun x H => proj1 (binary_normalize64_exact x H)) in H1 by smart_omega. do 2 rewrite (fun x H => proj2 (binary_normalize64_exact x H)) in H1 by smart_omega. - rewrite <- Z2R_plus, round_exact in H1 by smart_omega. - rewrite Rlt_bool_true in H1; + rewrite <- Z2R_plus, round_exact in H1 by smart_omega. + rewrite Rlt_bool_true in H1; replace (Int.unsigned x - 2147483648 + 2147483648) with (Int.unsigned x) in * by ring. apply B2R_inj. destruct (binary_normalize64_exact (Int.unsigned x)); [now smart_omega|]. @@ -616,7 +767,7 @@ Proof. now apply valid_rnd_round_mode. apply generic_format_FIX. exists (Float radix2 (cond_Zopp b (Zpos m)) 0). split; reflexivity. split; [reflexivity|]. - rewrite round_generic, Z2R_mult, Z2R_Zpower_pos, <- bpow_powerRZ; + rewrite round_generic, Z2R_mult, Z2R_Zpower_pos, <- bpow_powerRZ; [reflexivity|now apply valid_rnd_round_mode|apply generic_format_F2R; discriminate]. rewrite (inbetween_float_ZR_sign _ _ _ ((Zpos m) / Zpower_pos radix2 p) (new_location (Zpower_pos radix2 p) (Zpos m mod Zpower_pos radix2 p) loc_Exact)). @@ -685,7 +836,7 @@ Theorem intuoffloat_correct: end. Proof. intro; pose proof (Zoffloat_correct f); unfold intuoffloat; destruct (Zoffloat f). - pose proof (Zle_bool_spec 0 z); pose proof (Zle_bool_spec z Int.max_unsigned). + pose proof (Zle_bool_spec 0 z); pose proof (Zle_bool_spec z Int.max_unsigned). compute_this Int.max_unsigned; destruct H. inversion H0. inversion H1. rewrite <- (Int.unsigned_repr z) in H2 by smart_omega; split; assumption. @@ -713,14 +864,14 @@ Proof. simpl B2R; change 0%R with (Z2R 0); change (-1)%R with (Z2R (-1)); split; apply Z2R_lt; reflexivity. pose proof (Int.unsigned_range i). unfold round, scaled_mantissa, B2R, F2R, Fnum, Fexp in H0 |- *; simpl bpow in H0; do 2 rewrite Rmult_1_r in H0; - apply eq_Z2R in H0. + apply eq_Z2R in H0. split; apply Rnot_le_lt; intro. rewrite Ztrunc_ceil in H0; [apply Zceil_le in H3; change (-1)%R with (Z2R (-1)) in H3; rewrite Zceil_Z2R in H3; omega|]. eapply Rle_trans; [now apply H3|apply (Z2R_le (-1) 0); discriminate]. rewrite Ztrunc_floor in H0; [apply Zfloor_le in H3; rewrite Zfloor_Z2R in H3; now smart_omega|]. eapply Rle_trans; [|now apply H3]; apply (Z2R_le 0); discriminate. -Qed. +Qed. Theorem intuoffloat_intoffloat_1: forall x n, @@ -738,7 +889,7 @@ Proof. destruct H4; [rewrite H2 in H4; discriminate|]. apply intuoffloat_interval in H0; exfalso; destruct H0, H4. eapply Rlt_le_trans in H0; [|now apply H4]; apply (lt_Z2R (-1)) in H0; discriminate. - apply Rcompare_Lt_inv in H1; eapply Rle_lt_trans in H1; [|now apply H4]. + apply Rcompare_Lt_inv in H1; eapply Rle_lt_trans in H1; [|now apply H4]. unfold floatofintu in H1; rewrite (fun x H => proj1 (binary_normalize64_exact x H)) in H1; [apply lt_Z2R in H1; discriminate|split; reflexivity]. Qed. @@ -769,16 +920,16 @@ Proof. [rewrite H in H2; simpl B2R in H2; apply (eq_Z2R 0) in H2; discriminate|reflexivity]. reflexivity. rewrite H in H2; apply Rcompare_Gt_inv in H2; pose proof (intuoffloat_interval _ _ H1). - unfold sub, b64_minus. - exploit (Bminus_correct 53 1024 eq_refl eq_refl mode_NE x (floatofintu ox8000_0000)); [assumption|reflexivity|]; intro. + unfold sub, b64_minus. + exploit (Bminus_correct 53 1024 eq_refl eq_refl binop_pl mode_NE x (floatofintu ox8000_0000)); [assumption|reflexivity|]; intro. rewrite H, round_generic in H6. - match goal with [H6:if Rlt_bool ?x ?y then _ else _|-_] => + match goal with [H6:if Rlt_bool ?x ?y then _ else _|-_] => pose proof (Rlt_bool_spec x y); destruct (Rlt_bool x y) end. - destruct H6. + destruct H6 as [? []]. match goal with [|- _ ?y = _] => pose proof (intoffloat_correct y); destruct (intoffloat y) end. - destruct H9. + destruct H10. f_equal; rewrite <- (Int.repr_signed i); unfold Int.sub; f_equal; apply eq_Z2R. - rewrite Z2R_minus, H10, H4. + rewrite Z2R_minus, H11, H4. unfold round, scaled_mantissa, F2R, Fexp, Fnum, round_mode; simpl bpow; repeat rewrite Rmult_1_r; rewrite <- Z2R_minus; f_equal. rewrite (Ztrunc_floor (B2R _ _ x)), <- Zfloor_minus, <- Ztrunc_floor; @@ -786,11 +937,11 @@ Proof. left; eapply Rlt_trans; [|now apply H2]; apply (Z2R_lt 0); reflexivity. try (change (0 ?= 53) with Lt in H6,H8). (* for Coq 8.4 *) try (change (53 ?= 1024) with Lt in H6,H8). (* for Coq 8.4 *) - exfalso; simpl Zcompare in H6, H8; rewrite H6, H8 in H9. - destruct H9 as [|[]]; [discriminate|..]. - eapply Rle_trans in H9; [|apply Rle_0_minus; left; assumption]; apply (le_Z2R 0) in H9; apply H9; reflexivity. - eapply Rle_lt_trans in H9; [|apply Rplus_lt_compat_r; now apply (proj2 H5)]. - rewrite <- Z2R_opp, <- Z2R_plus in H9; apply lt_Z2R in H9; discriminate. + exfalso; simpl Zcompare in H6, H8; rewrite H6, H8 in H10. + destruct H10 as [|[]]; [discriminate|..]. + eapply Rle_trans in H10; [|apply Rle_0_minus; left; assumption]; apply (le_Z2R 0) in H10; apply H10; reflexivity. + eapply Rle_lt_trans in H10; [|apply Rplus_lt_compat_r; now apply (proj2 H5)]. + rewrite <- Z2R_opp, <- Z2R_plus in H10; apply lt_Z2R in H10; discriminate. exfalso; inversion H7; rewrite Rabs_right in H8. eapply Rle_lt_trans in H8. apply Rle_not_lt in H8; [assumption|apply (bpow_le _ 31); discriminate]. change (bpow radix2 31) with (Z2R(Zsucc Int.max_unsigned - Int.unsigned ox8000_0000)); rewrite Z2R_minus. @@ -816,22 +967,22 @@ Lemma split_bits_or: Proof. intros. transitivity (split_bits 52 11 (join_bits 52 11 false (Int.unsigned x) 1075)). - - f_equal. rewrite Int64.ofwords_add'. reflexivity. - - apply split_join_bits. + - f_equal. rewrite Int64.ofwords_add'. reflexivity. + - apply split_join_bits. compute; auto. - generalize (Int.unsigned_range x). - compute_this Int.modulus; compute_this (2^52); omega. + generalize (Int.unsigned_range x). + compute_this Int.modulus; compute_this (2^52); omega. compute_this (2^11); omega. Qed. Lemma from_words_value: forall x, - B2R _ _ (from_words ox4330_0000 x) = + B2R _ _ (from_words ox4330_0000 x) = (bpow radix2 52 + Z2R (Int.unsigned x))%R /\ is_finite _ _ (from_words ox4330_0000 x) = true. Proof. - intros; unfold from_words, double_of_bits, b64_of_bits, binary_float_of_bits. - rewrite B2R_FF2B; unfold is_finite; rewrite match_FF2B; + intros; unfold from_words, double_of_bits, b64_of_bits, binary_float_of_bits. + rewrite B2R_FF2B. rewrite is_finite_FF2B. unfold binary_float_of_bits_aux; rewrite split_bits_or; simpl; pose proof (Int.unsigned_range x). destruct (Int.unsigned x + Zpower_pos 2 52) eqn:?. exfalso; now smart_omega. @@ -848,37 +999,37 @@ Theorem floatofintu_from_words: sub (from_words ox4330_0000 x) (from_words ox4330_0000 Int.zero). Proof. intros; destruct (Int.eq_dec x Int.zero); [subst; vm_compute; reflexivity|]. - assert (Int.unsigned x <> 0). + assert (Int.unsigned x <> 0). intro; destruct n; rewrite <- (Int.repr_unsigned x), H; reflexivity. pose proof (Int.unsigned_range x). pose proof (binary_normalize64_exact (Int.unsigned x)). destruct H1; [smart_omega|]. unfold floatofintu, sub, b64_minus. - match goal with [|- _ = Bminus _ _ _ _ _ ?x ?y] => - pose proof (Bminus_correct 53 1024 eq_refl eq_refl mode_NE x y) end. + match goal with [|- _ = Bminus _ _ _ _ _ _ ?x ?y] => + pose proof (Bminus_correct 53 1024 eq_refl eq_refl binop_pl mode_NE x y) end. apply (fun f x y => f x y) in H3; try apply (fun x => proj2 (from_words_value x)). do 2 rewrite (fun x => proj1 (from_words_value x)) in H3. rewrite Int.unsigned_zero in H3. replace (bpow radix2 52 + Z2R (Int.unsigned x) - (bpow radix2 52 + Z2R 0))%R with (Z2R (Int.unsigned x)) in H3 by (simpl; ring). rewrite round_exact in H3 by smart_omega. - match goal with [H3:if Rlt_bool ?x ?y then _ else _ |- _] => - pose proof (Rlt_bool_spec x y); destruct (Rlt_bool x y) end; destruct H3. + match goal with [H3:if Rlt_bool ?x ?y then _ else _ |- _] => + pose proof (Rlt_bool_spec x y); destruct (Rlt_bool x y) end; destruct H3 as [? []]. try (change (53 ?= 1024) with Lt in H3,H5). (* for Coq 8.4 *) simpl Zcompare in *; apply B2R_inj; try match goal with [H':B2R _ _ ?f = _ , H'':is_finite _ _ ?f = true |- is_finite_strict _ _ ?f = true] => destruct f; [ simpl in H'; change 0%R with (Z2R 0) in H'; apply eq_Z2R in H'; now destruct (H (eq_sym H')) | discriminate H'' | discriminate H'' | reflexivity - ] + ] end. rewrite H3; assumption. - inversion H4; change (bpow radix2 1024) with (Z2R (radix2 ^ 1024)) in H6; rewrite <- Z2R_abs in H6. - apply le_Z2R in H6; pose proof (Zabs_spec (Int.unsigned x)); + inversion H4; change (bpow radix2 1024) with (Z2R (radix2 ^ 1024)) in H5; rewrite <- Z2R_abs in H5. + apply le_Z2R in H5; pose proof (Zabs_spec (Int.unsigned x)); exfalso; now smart_omega. Qed. Lemma ox8000_0000_signed_unsigned: - forall x, + forall x, Int.unsigned (Int.add x ox8000_0000) = Int.signed x + Int.half_modulus. Proof. intro; unfold Int.signed, Int.add; pose proof (Int.unsigned_range x). @@ -902,16 +1053,16 @@ Local Transparent Int.repr Int64.repr. pose proof (Int.signed_range x). pose proof (binary_normalize64_exact (Int.signed x)); destruct H1; [now smart_omega|]. unfold floatofint, sub, b64_minus. - match goal with [|- _ = Bminus _ _ _ _ _ ?x ?y] => - pose proof (Bminus_correct 53 1024 eq_refl eq_refl mode_NE x y) end. + match goal with [|- _ = Bminus _ _ _ _ _ _ ?x ?y] => + pose proof (Bminus_correct 53 1024 eq_refl eq_refl binop_pl mode_NE x y) end. apply (fun f x y => f x y) in H3; try apply (fun x => proj2 (from_words_value x)). do 2 rewrite (fun x => proj1 (from_words_value x)) in H3. replace (bpow radix2 52 + Z2R (Int.unsigned (Int.add x ox8000_0000)) - (bpow radix2 52 + Z2R (Int.unsigned ox8000_0000)))%R with (Z2R (Int.signed x)) in H3 by (rewrite ox8000_0000_signed_unsigned; rewrite Z2R_plus; simpl; ring). rewrite round_exact in H3 by smart_omega. - match goal with [H3:if Rlt_bool ?x ?y then _ else _ |- _] => - pose proof (Rlt_bool_spec x y); destruct (Rlt_bool x y) end; destruct H3. + match goal with [H3:if Rlt_bool ?x ?y then _ else _ |- _] => + pose proof (Rlt_bool_spec x y); destruct (Rlt_bool x y) end; destruct H3 as [? []]. try (change (0 ?= 53) with Lt in H3,H5). (* for Coq 8.4 *) try (change (53 ?= 1024) with Lt in H3,H5). (* for Coq 8.4 *) simpl Zcompare in *; apply B2R_inj; @@ -919,11 +1070,11 @@ Local Transparent Int.repr Int64.repr. destruct f; [ simpl in H'; change 0%R with (Z2R 0) in H'; apply eq_Z2R in H'; now destruct (H (eq_sym H')) | discriminate H'' | discriminate H'' | reflexivity - ] + ] end. rewrite H3; assumption. - inversion H4; unfold bpow in H6; rewrite <- Z2R_abs in H6; - apply le_Z2R in H6; pose proof (Zabs_spec (Int.signed x)); exfalso; now smart_omega. + inversion H4; unfold bpow in H5; rewrite <- Z2R_abs in H5; + apply le_Z2R in H5; pose proof (Zabs_spec (Int.signed x)); exfalso; now smart_omega. Qed. (** Conversions from 32-bit integers to single-precision floats can @@ -936,10 +1087,10 @@ Lemma is_finite_strict_ge_1: (1 <= Rabs (B2R _ _ f))%R -> is_finite_strict _ _ f = true. Proof. - intros. destruct f; auto. simpl in H0. + intros. destruct f; auto. simpl in H0. change 0%R with (Z2R 0) in H0. change 1%R with (Z2R 1) in H0. - rewrite <- Z2R_abs in H0. + rewrite <- Z2R_abs in H0. exploit le_Z2R; eauto. Qed. @@ -954,29 +1105,27 @@ Proof. subst n; reflexivity. exploit binary_normalize64_exact; eauto. intros [A B]. destruct (binary_normalize64 n 0 false) as [ | | | s m e] eqn:B64; simpl in *. -- assert (0 = n) by (apply eq_Z2R; auto). subst n. simpl in GT. omegaContradiction. +- assert (0 = n) by (apply eq_Z2R; auto). subst n. simpl in GT. omegaContradiction. - discriminate. - discriminate. - set (n1 := cond_Zopp s (Z.pos m)) in *. - generalize (binary_normalize32_correct n1 e s). + generalize (binary_normalize32_correct n1 e s). fold (binary_normalize32 n1 e s). intros C. - generalize (binary_normalize32_correct n 0 false). + generalize (binary_normalize32_correct n 0 false). fold (binary_normalize32 n 0 false). intros D. assert (A': @F2R radix2 {| Fnum := n; Fexp := 0 |} = Z2R n). - { - unfold F2R. apply Rmult_1_r. - } - rewrite A in C. rewrite A' in D. + { unfold F2R. apply Rmult_1_r. } + rewrite A in C. rewrite A' in D. destruct (Rlt_bool (Rabs (round radix2 (FLT_exp (3 - 128 - 24) 24) (round_mode mode_NE) (Z2R n))) (bpow radix2 128)). -+ destruct C as [C1 C2]; destruct D as [D1 D2]. ++ destruct C as [C1 [C2 _]]; destruct D as [D1 [D2 _]]. assert (1 <= Rabs (round radix2 (FLT_exp (3 - 128 - 24) 24) (round_mode mode_NE) (Z2R n)))%R. - { apply abs_round_ge_generic. + { apply abs_round_ge_generic. apply fexp_correct. red. omega. apply valid_rnd_round_mode. - apply generic_format_bpow with (e := 0). compute. congruence. + apply generic_format_bpow with (e := 0). compute. congruence. rewrite <- Z2R_abs. change 1%R with (Z2R 1). apply Z2R_le. omega. } apply B2R_inj. apply is_finite_strict_ge_1; auto. rewrite C1; auto. @@ -988,17 +1137,102 @@ Qed. Theorem singleofint_floatofint: forall n, singleofint n = singleoffloat (floatofint n). Proof. - intros. symmetry. apply single_float_of_int. + intros. symmetry. apply single_float_of_int. generalize (Int.signed_range n). smart_omega. Qed. Theorem singleofintu_floatofintu: forall n, singleofintu n = singleoffloat (floatofintu n). Proof. - intros. symmetry. apply single_float_of_int. + intros. symmetry. apply single_float_of_int. generalize (Int.unsigned_range n). smart_omega. Qed. +Theorem mul2_add: + forall f, add f f = mul f (floatofint (Int.repr 2%Z)). +Proof. + intros. unfold add, b64_plus, mul, b64_mult. + destruct (is_finite_strict _ _ f) eqn:EQFINST. + - assert (EQFIN:is_finite _ _ f = true) by (destruct f; simpl in *; congruence). + pose proof (Bplus_correct 53 1024 eq_refl eq_refl binop_pl mode_NE f f EQFIN EQFIN). + pose proof (Bmult_correct 53 1024 eq_refl eq_refl binop_pl mode_NE f + (floatofint (Int.repr 2%Z))). + rewrite <- double, Rmult_comm in H. + replace (B2R 53 1024 (floatofint (Int.repr 2))) with 2%R in H0 by (compute; field). + destruct Rlt_bool. + + destruct H0 as [? []], H as [? []]. + rewrite EQFIN in H1. + apply B2R_Bsign_inj; auto. + etransitivity. apply H. symmetry. apply H0. + etransitivity. apply H4. symmetry. etransitivity. apply H2. + destruct Bmult; try reflexivity; discriminate. + simpl. rewrite xorb_false_r. + erewrite <- Rmult_0_l, Rcompare_mult_r. + destruct f; try discriminate EQFINST. + simpl. unfold F2R. + erewrite <- Rmult_0_l, Rcompare_mult_r. + rewrite Rcompare_Z2R with (y:=0). + destruct b; reflexivity. + apply bpow_gt_0. + apply (Z2R_lt 0 2). omega. + + destruct H. + apply B2FF_inj. + etransitivity. apply H. + symmetry. etransitivity. apply H0. + f_equal. destruct Bsign; reflexivity. + - destruct f as [[]|[]| |]; try discriminate; try reflexivity. + simpl. erewrite binop_propagate1; reflexivity. +Qed. + +Program Definition pow2_float (b:bool) (e:Z) (H:-1023 < e < 1023) : float := + B754_finite _ _ b (nat_iter 52 xO xH) (e-52) _. +Next Obligation. + unfold Fappli_IEEE.bounded, canonic_mantissa. + rewrite andb_true_iff, Zle_bool_true by omega. split; auto. + apply Zeq_bool_true. unfold FLT_exp. simpl Z.of_nat. + apply Z.max_case_strong; omega. +Qed. + +Theorem mul_div_pow2: + forall b e f H H', + mul f (pow2_float b e H) = div f (pow2_float b (-e) H'). +Proof. + intros. unfold mul, b64_mult, div, b64_div. + pose proof (Bmult_correct 53 1024 eq_refl eq_refl binop_pl mode_NE f (pow2_float b e H)). + pose proof (Bdiv_correct 53 1024 eq_refl eq_refl binop_pl mode_NE f (pow2_float b (-e) H')). + lapply H1. clear H1. intro. + change (is_finite 53 1024 (pow2_float b e H)) with true in H0. + unfold Rdiv in H1. + replace (/ B2R 53 1024 (pow2_float b (-e) H'))%R + with (B2R 53 1024 (pow2_float b e H)) in H1. + destruct (is_finite _ _ f) eqn:EQFIN. + - destruct Rlt_bool. + + destruct H0 as [? []], H1 as [? []]. + apply B2R_Bsign_inj; auto. + etransitivity. apply H0. symmetry. apply H1. + etransitivity. apply H3. destruct Bmult; try discriminate H2; reflexivity. + symmetry. etransitivity. apply H5. destruct Bdiv; try discriminate H4; reflexivity. + reflexivity. + + apply B2FF_inj. + etransitivity. apply H0. symmetry. etransitivity. apply H1. + reflexivity. + - destruct f; try discriminate EQFIN. reflexivity. + simpl. erewrite binop_propagate2. reflexivity. + reflexivity. reflexivity. reflexivity. + - simpl. + assert ((4503599627370496 * bpow radix2 (e - 52))%R = + (/ (4503599627370496 * bpow radix2 (- e - 52)))%R). + { etransitivity. symmetry. apply (bpow_plus radix2 52). + symmetry. etransitivity. apply f_equal. symmetry. apply (bpow_plus radix2 52). + rewrite <- bpow_opp. f_equal. ring. } + destruct b. unfold cond_Zopp. + rewrite !F2R_Zopp, <- Ropp_inv_permute. f_equal. auto. + intro. apply F2R_eq_0_reg in H3. omega. + apply H2. + - simpl. intro. apply F2R_eq_0_reg in H2. + destruct b; simpl in H2; omega. +Qed. + Global Opaque zero eq_dec neg abs singleoffloat intoffloat intuoffloat floatofint floatofintu add sub mul div cmp bits_of_double double_of_bits bits_of_single single_of_bits from_words. diff --git a/powerpc/Nan.v b/powerpc/Nan.v new file mode 100644 index 0000000..77a752e --- /dev/null +++ b/powerpc/Nan.v @@ -0,0 +1,28 @@ +Require Import Fappli_IEEE. +Require Import Fappli_IEEE_bits. +Require Import Floats. +Require Import ZArith. +Require Import Integers. + +(* Needed to break a circular reference after extraction *) +Definition transform_quiet_pl := + Eval unfold Float.transform_quiet_pl in Float.transform_quiet_pl. + +Program Definition default_pl : bool * nan_pl 53 := (false, nat_iter 51 xO xH). + +Definition binop_pl (pl1 pl2:binary64) : bool*nan_pl 53 := + match pl1, pl2 with + | B754_nan s1 pl1, _ => (s1, transform_quiet_pl pl1) + | _, B754_nan s2 pl2 => (s2, transform_quiet_pl pl2) + | _, _ => default_pl + end. + +Theorem binop_propagate1: Float.binop_propagate1_prop binop_pl. +Proof. + repeat intro. destruct f1, f2; try discriminate; reflexivity. +Qed. + +Theorem binop_propagate2: Float.binop_propagate2_prop binop_pl. +Proof. + repeat intro. destruct f1, f2, f3; try discriminate; simpl; reflexivity. +Qed. -- cgit v1.2.3