1# This file is part of crlibm, the correctly rounded mathematical library,
2# which has been developed by the Arénaire project at École normale supérieure
3# de Lyon.
4#
5# Copyright (C) 2004-2011 David Defour, Catherine Daramy-Loirat,
6# Florent de Dinechin, Matthieu Gallet, Nicolas Gast, Christoph Quirin Lauter,
7# and Jean-Michel Muller
8#
9# This library is free software; you can redistribute it and/or
10# modify it under the terms of the GNU Lesser General Public
11# License as published by the Free Software Foundation; either
12# version 2.1 of the License, or (at your option) any later version.
13#
14# This library is distributed in the hope that it will be useful,
15# but WITHOUT ANY WARRANTY; without even the implied warranty of
16# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
17# Lesser General Public License for more details.
18#
19# You should have received a copy of the GNU Lesser General Public
20# License along with this library; if not, write to the Free Software
21# Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
22
23accPolyC14 = <float64ne> (_accPolyC14);
24accPolyC13 = <float64ne> (_accPolyC13);
25accPolyC12 = <float64ne> (_accPolyC12);
26accPolyC11 = <float64ne> (_accPolyC11);
27accPolyC10 = <float64ne> (_accPolyC10);
28
29accPolyC9h = <float64ne> (_accPolyC9h);
30accPolyC9l = <float64ne> (_accPolyC9l);
31AccPolyC9hl = accPolyC9h + accPolyC9l; #MAPLE
32accPolyC8h = <float64ne> (_accPolyC8h);
33accPolyC8l = <float64ne> (_accPolyC8l);
34AccPolyC8hl = accPolyC8h + accPolyC8l; #MAPLE
35accPolyC7h = <float64ne> (_accPolyC7h);
36accPolyC7l = <float64ne> (_accPolyC7l);
37AccPolyC7hl = accPolyC7h + accPolyC7l; #MAPLE
38accPolyC6h = <float64ne> (_accPolyC6h);
39accPolyC6l = <float64ne> (_accPolyC6l);
40AccPolyC6hl = accPolyC6h + accPolyC6l; #MAPLE
41accPolyC5h = <float64ne> (_accPolyC5h);
42accPolyC5l = <float64ne> (_accPolyC5l);
43AccPolyC5hl = accPolyC5h + accPolyC5l; #MAPLE
44accPolyC4h = <float64ne> (_accPolyC4h);
45accPolyC4l = <float64ne> (_accPolyC4l);
46AccPolyC4hl = accPolyC4h + accPolyC4l; #MAPLE
47accPolyC3h = <float64ne> (_accPolyC3h);
48accPolyC3l = <float64ne> (_accPolyC3l);
49AccPolyC3hl = accPolyC3h + accPolyC3l; #MAPLE
50
51E = 0; #MAPLE
52
53zh = <float64ne> (Z);
54zl = Z - zh; #MAPLE
55
56highPoly <float64ne> = accPolyC10 + zh * (accPolyC11 + zh * (accPolyC12 + zh * (accPolyC13 + zh * accPolyC14)));
57
58T1hl = zh * highPoly; #MAPLE
59
60T2 = AccPolyC9hl + T1hl; #MAPLE
61T3 = Z * T2hl; #MAPLE
62T4 = AccPolyC8hl + T3hl; #MAPLE
63T5 = Z * T4hl; #MAPLE
64T6 = AccPolyC7hl + T5hl; #MAPLE
65T7 = Z * T6hl; #MAPLE
66T8 = AccPolyC6hl + T7hl; #MAPLE
67T9 = Z * T8hl; #MAPLE
68T10 = AccPolyC5hl + T9hl; #MAPLE
69T11 = Z * T10hl; #MAPLE
70T12 = AccPolyC4hl + T11hl; #MAPLE
71T13 = Z * T12hl; #MAPLE
72T14 = AccPolyC3hl + T13hl; #MAPLE
73
74
75ZSquare = Z * Z; #MAPLE
76ZCube = Z * ZSquarehml; #MAPLE
77HigherPolyMultZ = T14hl * ZCubehml; #MAPLE
78ZSquareHalfhml = -0.5 * ZSquarehml; #MAPLE
79PolyWithSquare = ZSquareHalfhml + HigherPolyMultZhml; #MAPLE
80Poly = Z + PolyWithSquarehml; #MAPLE
81
82#We can simplify the proof in this case since we know that adding a triple double which is
83#equal to 0 exactly is exact.
84
85Loghml = Polyhml; #MAPLE
86
87
88#Mathematical definition of the logarithm
89
90MHighPoly = accPolyC10 + Z * (accPolyC11 + Z * (accPolyC12 + Z * (accPolyC13 + Z * accPolyC14))); #MAPLE
91MT1 = Z * MHighPoly; #MAPLE
92MT2 = AccPolyC9hl + MT1; #MAPLE
93MT3 = Z * MT2; #MAPLE
94MT4 = AccPolyC8hl + MT3; #MAPLE
95MT5 = Z * MT4; #MAPLE
96MT6 = AccPolyC7hl + MT5; #MAPLE
97MT7 = Z * MT6; #MAPLE
98MT8 = AccPolyC6hl + MT7; #MAPLE
99MT9 = Z * MT8; #MAPLE
100MT10 = AccPolyC5hl + MT9; #MAPLE
101MT11 = Z * MT10; #MAPLE
102MT12 = AccPolyC4hl + MT11; #MAPLE
103MT13 = Z * MT12; #MAPLE
104MT14 = AccPolyC3hl + MT13; #MAPLE
105MZSquare = Z * Z; #MAPLE
106MZCube = Z * MZSquare; #MAPLE
107MHigherPolyMultZ = MT14 * MZCube; #MAPLE
108MZSquareHalf = -0.5 * MZSquare; #MAPLE
109MPolyWithSquare = MZSquareHalf + MHigherPolyMultZ; #MAPLE
110MPoly = Z + MPolyWithSquare; #MAPLE
111
112#We apply the same simplification here
113
114MLog = MLog1pZ; #MAPLE
115
116
117#Useful additional definitions
118
119epsilon1 = (highPoly - MHighPoly) / MHighPoly; #MAPLE
120epsilon2 = (T1hl - MT1) / MT1; #MAPLE
121epsilon3 = (T2hl - MT2) / MT2; #MAPLE
122epsilon4 = (T3hl - MT3) / MT3; #MAPLE
123epsilon5 = (T4hl - MT4) / MT4; #MAPLE
124epsilon6 = (T5hl - MT5) / MT5; #MAPLE
125epsilon7 = (T6hl - MT6) / MT6; #MAPLE
126epsilon8 = (T7hl - MT7) / MT7; #MAPLE
127epsilon9 = (T8hl - MT8) / MT8; #MAPLE
128epsilon10 = (T9hl - MT9) / MT9; #MAPLE
129epsilon11 = (T10hl - MT10) / MT10; #MAPLE
130epsilon12 = (T11hl - MT11) / MT11; #MAPLE
131epsilon13 = (T12hl - MT12) / MT12; #MAPLE
132epsilon14 = (T13hl - MT13) / MT13; #MAPLE
133epsilon15 = (T14hl - MT14) / MT14; #MAPLE
134
135epsilon16 = (ZCubehml - MZCube) / MZCube; #MAPLE
136epsilon17 = (HigherPolyMultZhml - MHigherPolyMultZ) / MHigherPolyMultZ; #MAPLE
137epsilon18 = (ZSquareHalfhml - MZSquareHalf) / MZSquareHalf; #MAPLE
138epsilon19 = (PolyWithSquarehml - MPolyWithSquare) / MPolyWithSquare; #MAPLE
139epsilon20 = (Polyhml - MLog1pZ) / MLog1pZ; #MAPLE
140
141epsilon21 = (PolyWithSquare - MPolyWithSquare) / MPolyWithSquare; #MAPLE
142epsilon22 = (Polyhml - MPoly) / MPoly; #MAPLE
143epsilon23 = (Poly - MPoly) / MPoly; #MAPLE
144
145aux1 = -0.5 * Z + MZSquare * MT14; #MAPLE
146
147
148#End additional definitions
149
150{
151(
152(T2hl - T2) / T2 in [-1b-103,1b-103]
153/\ (T3hl - T3) / T3 in [-1b-102,1b-102]
154/\ (T4hl - T4) / T4 in [-1b-103,1b-103]
155/\ (T5hl - T5) / T5 in [-1b-102,1b-102]
156/\ (T6hl - T6) / T6 in [-1b-103,1b-103]
157/\ (T7hl - T7) / T7 in [-1b-102,1b-102]
158/\ (T8hl - T8) / T8 in [-1b-103,1b-103]
159/\ (T9hl - T9) / T9 in [-1b-102,1b-102]
160/\ (T10hl - T10) / T10 in [-1b-103,1b-103]
161/\ (T11hl - T11) / T11 in [-1b-102,1b-102]
162/\ (T12hl - T12) / T12 in [-1b-103,1b-103]
163/\ (T13hl - T13) / T13 in [-1b-102,1b-102]
164/\ (T14hl - T14) / T14 in [-1b-103,1b-103]
165/\ (ZSquarehml - ZSquare) / ZSquare in [-1b-149,1b-149]
166/\ (ZCubehml - ZCube) / ZCube in [-1b-144,1b-144]
167/\ (HigherPolyMultZhml - HigherPolyMultZ) / HigherPolyMultZ in [-1b-141,1b-141]
168/\ (PolyWithSquarehml - PolyWithSquare) / PolyWithSquare in [-1b-137,1b-137]
169/\ (Polyhml - Poly) / Poly in [-1b-134,1b-134]
170/\ (MPoly - MLog1pZ) / MLog1pZ in [-_epsilonApproxAccurate,_epsilonApproxAccurate]
171/\ Z in [1b-900,_zmax]
172/\ ((logh + logm + logl) - Loghml) / Loghml in [-1b-159,1b-159]
173->
174((logh + logm + logl) - MLog) / MLog in [-5735b-132,5735b-132]
175)
176/\
177(
178(T2hl - T2) / T2 in [-1b-103,1b-103]
179/\ (T3hl - T3) / T3 in [-1b-102,1b-102]
180/\ (T4hl - T4) / T4 in [-1b-103,1b-103]
181/\ (T5hl - T5) / T5 in [-1b-102,1b-102]
182/\ (T6hl - T6) / T6 in [-1b-103,1b-103]
183/\ (T7hl - T7) / T7 in [-1b-102,1b-102]
184/\ (T8hl - T8) / T8 in [-1b-103,1b-103]
185/\ (T9hl - T9) / T9 in [-1b-102,1b-102]
186/\ (T10hl - T10) / T10 in [-1b-103,1b-103]
187/\ (T11hl - T11) / T11 in [-1b-102,1b-102]
188/\ (T12hl - T12) / T12 in [-1b-103,1b-103]
189/\ (T13hl - T13) / T13 in [-1b-102,1b-102]
190/\ (T14hl - T14) / T14 in [-1b-103,1b-103]
191/\ (ZSquarehml - ZSquare) / ZSquare in [-1b-149,1b-149]
192/\ (ZCubehml - ZCube) / ZCube in [-1b-144,1b-144]
193/\ (HigherPolyMultZhml - HigherPolyMultZ) / HigherPolyMultZ in [-1b-141,1b-141]
194/\ (PolyWithSquarehml - PolyWithSquare) / PolyWithSquare in [-1b-137,1b-137]
195/\ (Polyhml - Poly) / Poly in [-1b-134,1b-134]
196/\ (MPoly - MLog1pZ) / MLog1pZ in [-_epsilonApproxAccurate,_epsilonApproxAccurate]
197/\ Z in [_zmin,-1b-900]
198/\ ((logh + logm + logl) - Loghml) / Loghml in [-1b-159,1b-159]
199->
200((logh + logm + logl) - MLog) / MLog in [-5735b-132,5735b-132]
201)
202}
203
204((logh + logm + logl) - MLog) / MLog -> ((Loghml - MLog) / MLog) + ((((logh + logm + logl) - Loghml) / Loghml) * (((Loghml - MLog) / MLog) + 1));
205
206T2hl -> (T2 * ((T2hl - T2) / T2)) + T2;
207T3hl -> (T3 * ((T3hl - T3) / T3)) + T3;
208T4hl -> (T4 * ((T4hl - T4) / T4)) + T4;
209T5hl -> (T5 * ((T5hl - T5) / T5)) + T5;
210T6hl -> (T6 * ((T6hl - T6) / T6)) + T6;
211T7hl -> (T7 * ((T7hl - T7) / T7)) + T7;
212T8hl -> (T8 * ((T8hl - T8) / T8)) + T8;
213T9hl -> (T9 * ((T9hl - T9) / T9)) + T9;
214T10hl -> (T10 * ((T10hl - T10) / T10)) + T10;
215T11hl -> (T11 * ((T11hl - T11) / T11)) + T11;
216T12hl -> (T12 * ((T12hl - T12) / T12)) + T12;
217T13hl -> (T13 * ((T13hl - T13) / T13)) + T13;
218T14hl -> (T14 * ((T14hl - T14) / T14)) + T14;
219
220
221ZSquarehml -> (ZSquare * ((ZSquarehml - ZSquare) / ZSquare)) + ZSquare;
222ZCubehml -> (ZCube * ((ZCubehml - ZCube) / ZCube)) + ZCube;
223HigherPolyMultZhml -> (HigherPolyMultZ * ((HigherPolyMultZhml - HigherPolyMultZ) / HigherPolyMultZ)) + HigherPolyMultZ;
224PolyWithSquarehml -> (PolyWithSquare * ((PolyWithSquarehml - PolyWithSquare) / PolyWithSquare)) + PolyWithSquare;
225Polyhml -> (Poly * ((Polyhml - Poly) / Poly)) + Poly;
226
227
228epsilon2 -> epsilon1 + (((zh - Z) / Z) * (epsilon1 + 1));
229epsilon3 -> ((epsilon2 * MT1) / (AccPolyC9hl + MT1)) + (((AccPolyC9hl + T1hl) / (AccPolyC9hl + MT1)) * ((T2hl - T2) / T2));
230epsilon4 -> epsilon3 + (((T3hl - T3) / T3) * (epsilon3 + 1));
231epsilon5 -> ((epsilon4 * MT3) / (AccPolyC8hl + MT3)) + (((AccPolyC8hl + T3hl) / (AccPolyC8hl + MT3)) * ((T4hl - T4) / T4));
232epsilon6 -> epsilon5 + (((T5hl - T5) / T5) * (epsilon5 + 1));
233epsilon7 -> ((epsilon6 * MT5) / (AccPolyC7hl + MT5)) + (((AccPolyC7hl + T5hl) / (AccPolyC7hl + MT5)) * ((T6hl - T6) / T6));
234epsilon8 -> epsilon7 + (((T7hl - T7) / T7) * (epsilon7 + 1));
235epsilon9 -> ((epsilon8 * MT7) / (AccPolyC6hl + MT7)) + (((AccPolyC6hl + T7hl) / (AccPolyC6hl + MT7)) * ((T8hl - T8) / T8));
236epsilon10 -> epsilon9 + (((T9hl - T9) / T9) * (epsilon9 + 1));
237epsilon11 -> ((epsilon10 * MT9) / (AccPolyC5hl + MT9)) + (((AccPolyC5hl + T9hl) / (AccPolyC5hl + MT9)) * ((T10hl - T10) / T10));
238epsilon12 -> epsilon11 + (((T11hl - T11) / T11) * (epsilon11 + 1));
239epsilon13 -> ((epsilon12 * MT11) / (AccPolyC4hl + MT11)) + (((AccPolyC4hl + T11hl) / (AccPolyC4hl + MT11)) * ((T12hl - T12) / T12));
240epsilon14 -> epsilon13 + (((T13hl - T13) / T13) * (epsilon13 + 1));
241epsilon15 -> ((epsilon14 * MT13) / (AccPolyC3hl + MT13)) + (((AccPolyC3hl + T13hl) / (AccPolyC3hl + MT13)) * ((T14hl - T14) / T14));
242epsilon16 -> ((ZSquarehml - MZSquare) / MZSquare) + (((ZCubehml - ZCube) / ZCube) * (((ZSquarehml - MZSquare) / MZSquare) + 1));
243epsilon17 -> epsilon15 + epsilon16 + epsilon15 * epsilon16 +
244((HigherPolyMultZhml - HigherPolyMultZ) / HigherPolyMultZ) * (1 + epsilon15 + epsilon16 + epsilon15 * epsilon16);
245epsilon18 -> (ZSquarehml - MZSquare) / MZSquare;
246epsilon19 -> epsilon21 + (1 + epsilon21) * ((PolyWithSquarehml - PolyWithSquare) / PolyWithSquare);
247epsilon20 -> (((Polyhml - MPoly) / MPoly) + ((MPoly - MLog1pZ) / MLog1pZ)) + (((Polyhml - MPoly) / MPoly) * ((MPoly - MLog1pZ) / MLog1pZ));
248epsilon21 -> ((-0.5 * epsilon18) + (Z * MT14 * epsilon17)) / (-0.5 + (Z * MT14));
249epsilon22 -> epsilon23 + (((Polyhml - Poly) / Poly) * (1+ epsilon23));
250epsilon23 -> epsilon19 * (aux1 / (1 + aux1));
251