1
2
3\environment lowlevelstyle
4
5\startdocument
6 [title=accuracy,
7 color=darkgray]
8
9\startsectionlevel[title=Introduction]
10
11When you look at \TEX\ and \METAPOST\ output the accuracy of the rendering stands
12out, unless of course you do a sloppy job on design and interfere badly with the
13system. Much has to do with the fact that calculations are very precise,
14especially given the time when \TEX\ was written. Because \TEX\ doesnt rely on
15(at that time nonportable) floating point calculations, it does all with 32
16bit integers, except in the backend where glue calculations are used for
17finalizing the glue values. It all changed a bit when we added \LUA\ because
18there we mix integers and doubles but in practice it works out okay.
19
20When looking at floating point (and posits) one can end up in discussions about
21which one is better, what the flaws fo each are, etc. Here were only interested
22in the fact that posits are more accurate in the ranges where \TEX\ and
23\METAPOST\ operate, as well as the fact that we only have 32 bits for floats in
24\TEX, unless we patch more heavily. So, it is also very much about storage.
25
26When you work with dimensions like points, they get converted to an integer
27number (the \type {sp} unit) and from that its just integer calculations. The
28maximum dimension is \the\maxdimen, which already shows a rounding issue. Of
29course when one goes precise for sure there is some loss, but on the average
30were okay. So, in the next example the two last rows are equivalent:
31
32\starttabulate[Trrr]
33\NC .1pt \NC \the\dimexpr.1pt \relax \NC \number\dimexpr.1pt \relax sp \NC \NR
34\NC .2pt \NC \the\dimexpr.2pt \relax \NC \number\dimexpr.2pt \relax sp \NC \NR
35\NC .3pt \NC \the\dimexpr.3pt \relax \NC \number\dimexpr.3pt \relax sp \NC \NR
36\NC .1pt .2pt \NC \the\dimexpr.1pt.2pt\relax \NC \number\dimexpr.1pt.2pt\relax sp \NC \NR
37\stoptabulate
38
39When were at the \LUA\ end things are different, there numbers are mapped onto
4064 bit floating point variables (doubles) and not all numbers map well. This is
41what we get when we work with doubles in \LUA:
42
43\starttabulate[Trr]
44\NC .1 \NC \luaexpr{.1 } \NC \NR
45\NC .2 \NC \luaexpr{.2 } \NC \NR
46\NC .3 \NC \luaexpr{.3 } \NC \NR
47\NC .1 .2 \NC \luaexpr{.1+.2} \NC \NR
48\stoptabulate
49
50The serialization looks as if all is okay but when we test for equality there
51is a problem:
52
53\starttabulate[Trl]
54\NC .3 == .3 \NC \luaexpr{tostring( .3 == .3)} \NC \NR
55\NC .1 .2 == .3 \NC \luaexpr{tostring(.1 + .2 == .3)} \NC \NR
56\stoptabulate
57
58This means that a test like this can give false positives or negatives unless one
59tests the difference against the accuracy (in \METAPOST\ we have the {eps}
60variable for that). In \TEX\ clipping of the decimal fraction influences equality.
61
62\starttabulate[Trl]
63\NC \type{\iflua { .3 == .3 } Y\else N\fi} \NC \iflua{ .3 == .3} equal\else different\fi \NC \NR
64\NC \type{\iflua { .1 .2 == .3 } Y\else N\fi} \NC \iflua{.1 .2 == .3} equal\else different\fi \NC \NR
65\stoptabulate
66
67The serialization above misguides us because the number of digits displayed is
68limited. Actually, when we would compare serialized strings the equality holds,
69definitely within the accuracy of \TEX. But here is reality:
70
71\startluacode
72 local a = 0.1
73 local b = 0.2
74 local c = 0.3
75 local d = a + b
76 local NC, NR = context.NC, context.NR
77 local function show(f)
78 context.NC() context(context.escaped(f))
79 context.NC() context(f,c)
80 context.NC() context(f,d)
81 context.NC() context.NR()
82 end
83 context.starttabulate { "|T|l|l|" }
84 context.NC()
85 context.NC() context(".3")
86 context.NC() context(".1 + .2")
87 context.NC() context.NR()
88
89 show("%0.10g",c)
90 show("%0.17g",c)
91 show("%0.20g",c)
92 show("%0.25g",c)
93 context.stoptabulate()
94\stopluacode
95
96The above examples use $0.1$, $0.2$ and $0.3$ and on a 32 bit float that actually
97works out okay, but \LUAMETATEX\ is 64 bit. Is this really important in practice?
98There are indeed cases where we are bitten by this. At the \LUA\ end we seldom
99test for equality on calculated values but it might impact check for less or
100greater then. At the \TEX\ end there are a few cases where we have issues but
101these also relate to the limited precision. It is not uncommon to loose a few
102scaled points so that has to be taken into account then. So how can we deal with
103this? In the next section(s) an alternative approach is discussed. It is not so
104much the solution for all problems but who knows.
105
106\stopsectionlevel
107
108\startsectionlevel[title=Posits]
109
110
111
112The next table shows the same as what we started with but with a different
113serialization.
114
115\starttabulate[Trr]
116\NC .1 \NC \positunum{.1 } \NC \NR
117\NC .2 \NC \positunum{.2 } \NC \NR
118\NC .3 \NC \positunum{.3 } \NC \NR
119\NC .1 .2 \NC \positunum{.1 .2} \NC \NR
120\stoptabulate
121
122And here we get equality in both cases:
123
124\starttabulate[Trl]
125\NC .3 == .3 \NC \positunum{ .3 == .3} \NC \NR
126\NC .1 .2 == .3 \NC \positunum{.1 .2 == .3} \NC \NR
127\stoptabulate
128
129The next table shows what we actually are dealing with. The \type {\if}test is
130not a primitive but provided by \CONTEXT.
131
132\starttabulate[Trl]
133\NC \type{\ifpositunum { .3 == .3 } Y\else N\fi} \NC \ifpositunum{ .3 == .3} equal\else different\fi \NC \NR
134\NC \type{\ifpositunum { .1 .2 == .3 } Y\else N\fi} \NC \ifpositunum{.1 .2 == .3} equal\else different\fi \NC \NR
135\stoptabulate
136
137And what happens when we do more complex calculations:
138
139\starttabulate[Trl]
140\NC \type {math .sin(0.1 0.2) == math .sin(0.3)} \NC \luaexpr{tostring(math.sin(0.1 + 0.2) == math.sin(0.3))} \NC \NR
141\NC \type {posit.sin(0.1 0.2) == posit.sin(0.3)} \NC \positunum{sin(0.1 0.2) == sin(0.3)} \NC \NR
142\stoptabulate
143
144Of course other numbers might work out differently! I just took the simple tests
145that came to mind.
146
147So what are these posits? Here its enough to know that they are a different way
148to store numbers with fractions. They still can loose precision but a bit less on
149smaller values and often we have relative small values in \TEX. Here are some links:
150
151\starttyping
152https:www.johngustafson.netpdfsBeatingFloatingPoint.pdf
153https:posithub.orgconga2019docs141130FlorentDeDinechin.pdf
154\stoptyping
155
156There are better explanations out there than I can provide (if at all). When I
157first read about these unums (a review of the 2015 book \quotation {The End of
158Error Unum Computing}) I was intrigued and when in 2023 I read something about it
159in relation to RISCV I decided to just add this playground for the users. After
160all we also have decimal support. And interval based solutions might actually be
161good for \METAPOST, so that is why we have it as extra number model. There we
162need to keep in mind that \METAPOST\ in non scaled models also apply some of the
163range checking and clipping that happens in scaled (these magick 4096 tricks).
164
165For now it is enough to know that its an alternative for floats that {\em could}
166work better in some cases but not all. The presentation mentioned above gives
167some examples of physics constants where 32 posits are not good enough for
168encoding the extremely large or small constants, but for $\pi$ its all fine.
169\footnote {Are 64 bit posits actually being worked on in softposit? There are
170some commented sections. We also need to patch some unions to make it compile as
171C.} In double mode we actually have quite good precision compared to 32 bit
172posits but with 32 bit floats we gain some. Very small numbers and very large
173numbers are less precise, but around 1 we gain: the next value after 1 is
1741.0000001 for a float and 1.000000008 for a posit (both 32 bit). So, currently
175for \METAPOST\ there is no real gain but if wed add posits to \TEX\ we could
176gain some because there a halfword (used for storing data) is 32 bit.
177
178But how about \TEX ? Per April 2023 the \LUAMETATEX\ engine has native support
179for floats (this in addition to \LUA\ based floats that we already had in
180\CONTEXT). How that works can be demonstrated with some examples. The float
181related commands are similar to those for numbers and dimensions: \typ
182{\floatdef}, \typ {\float}, \typ {\floatexpr}, \typ {\iffloat}, \typ
183{\ifzerofloat} and \typ {\ifintervalfloat}. That means that we also have them as
184registers. The \typ {\positdef} primitive is similar to \typ {\dimensiondef}.
185When a float (posit) is seen in a dimension context it will be interpreted as
186points, and in an integer context it will be a rounded number. As with other
187registers we have a \typ {\newfloat} macro. The \typ {\advance}, \typ
188{\multiply} and \typ {\divide} primitives accept floats.
189
190\startbuffer[reset]
191\scratchdimen=1.23456pt
192\scratchfloat=1.23456
193\stopbuffer
194
195\typebuffer[reset] \getbuffer[reset]
196
197We now use these two variables in an example:
198
199\startbuffer[dimensions]
200\setbox0\hbox to \scratchdimen {x}\the\wd0
201\scratchdimen \dimexpr \scratchdimen * 2\relax
202\setbox0\hbox to \scratchdimen {x}\the\wd0
203\advance \scratchdimen \scratchdimen
204\setbox0\hbox to \scratchdimen {x}\the\wd0
205\multiply\scratchdimen by 2
206\setbox0\hbox to \scratchdimen {x}\the\wd0
207\stopbuffer
208
209\typebuffer[dimensions] \startlines\darkblue\tttf\getbuffer[reset,dimensions]\stoplines
210
211When we use floats we get this:
212
213\startbuffer[floats]
214\setbox0\hbox to \scratchfloat {x}\the\wd0
215\scratchfloat \floatexpr \scratchfloat * 2\relax
216\setbox0\hbox to \scratchfloat {x}\the\wd0
217\advance \scratchfloat \scratchfloat
218\setbox0\hbox to \scratchfloat {x}\the\wd0
219\multiply\scratchfloat by 2
220\setbox0\hbox to \scratchfloat {x}\the\wd0
221\stopbuffer
222
223\typebuffer[floats] \startlines\darkblue\tttf\getbuffer[reset,floats]\stoplines
224
225So which approach is more accurate? At first sight you might think that the
226dimensions are better because in the last two lines they indeed duplicate.
227However, the next example shows that with dimensions we lost some between steps.
228
229\startbuffer[noboxes]
230 \the\scratchfloat
231\scratchfloat \floatexpr \scratchfloat * 2\relax \the\scratchfloat
232\advance \scratchfloat \scratchfloat \the\scratchfloat
233\multiply\scratchfloat by 2 \the\scratchfloat
234\stopbuffer
235
236\typebuffer[noboxes] \startlines\darkblue\tttf\getbuffer[reset,noboxes]\stoplines
237
238One problem with accuracy is that it can build up. So when one eventually does
239some comparison the expectations can be wrong.
240
241\startbuffer
242\dimen0=1.2345pt
243\dimen2=1.2345pt
244
245\ifdim \dimen0=\dimen2 S\else D\fi \space 0sp: [dim]
246\ifintervaldim0sp\dimen0 \dimen2 O\else D\fi \space 0sp: [0sp]
247
248\advance\dimen2 1sp
249
250\ifdim \dimen0=\dimen2 S\else D\fi \space 1sp: [dim]
251\ifintervaldim 1sp \dimen0 \dimen2 O\else D\fi \space 1sp: [1sp]
252\ifintervaldim 1sp \dimen2 \dimen0 O\else D\fi \space 1sp: [1sp]
253\ifintervaldim 2sp \dimen0 \dimen2 O\else D\fi \space 1sp: [2sp]
254\ifintervaldim 2sp \dimen2 \dimen0 O\else D\fi \space 1sp: [2sp]
255
256\advance\dimen2 1sp
257
258\ifintervaldim 1sp \dimen0\dimen2 O\else D\fi \space 2sp: [1sp]
259\ifintervaldim 1sp \dimen2\dimen0 O\else D\fi \space 2sp: [1sp]
260\ifintervaldim 5sp \dimen0\dimen2 O\else D\fi \space 2sp: [5sp]
261\ifintervaldim 5sp \dimen2\dimen0 O\else D\fi \space 2sp: [5sp]
262\stopbuffer
263
264\typebuffer
265
266Here we show a test for overlap in values, the same can be done with integer
267numbers (counts) and floats. This interval checking is an experiment and well
268see it if gets used.
269
270\startpacked\darkblue \tttf \getbuffer \stoppacked
271
272There are also \typ {\ifintervalfloat} and \typ{\ifintervalnum}. Because I have
273worked around these few scaled point rounding issues for decades, it might
274actually take some time before we see the interval tests being used in \CONTEXT.
275After all, there is no reason to touch somewhat tricky mechanism without reason
276(read: users complaining).
277
278To come back to posits, just to be clear, we use 32 bit posits and not 32 bit
279floats, which we could have but that way we gain some accuracy because less bits
280are used by default for the exponential.
281
282In \CONTEXT\ we also provide a bunch of pseudo primitives. These take one float:
283\type {\pfsin}, \type {\pfcos}, \type {\pftan}, \type {\pfasin}, \type {\pfacos},
284\type {\pfatan}, \type {\pfsinh}, \type {\pfcosh}, \type {\pftanh}, \type
285{\pfasinh}, \type {\pfacosh}, \type {\pfatanh}, \type {\pfsqrt}, \type {\pflog},
286\type {\pfexp}, \type {\pfceil}, \type {\pffloor}, \type {\pfround}, \type
287{\pfabs}, \type {\pfrad} and \type {\pfdeg}, whiel these expect two floats: \type
288{\pfatantwo}, \type {\pfpow}, \type {\pfmod} and \type {\pfrem}.
289
290
291
292\stopsectionlevel
293
294\startsectionlevel[title=\METAPOST]
295
296In addition to the instances \typ {metafun} (double in \LMTX), \typ {scaledfun},
297\typ {doublefun}, \typ {decimalfun} we now also have \typ {positfun}. Because we
298currently use 32 bit posits in the new number system there is no real gain over
299the already present 64 bit doubles. When 64 bit posits show up we might move on
300to that.
301
302\stopsectionlevel
303
304\startsectionlevel[title=\LUA]
305
306We support posits in \LUA\ too. Here we need to create a posit user data
307object. The usual metatable magick kicks in:
308
309\starttyping
310local p = posit.new(123.456)
311local q = posit.new(789.123)
312local r = p q
313\stoptyping
314
315Here we just mention what is currently interface. The management functions are:
316\typ {new}, \type {copy}, \typ {tostring}, \typ {tonumber}, \typ {integer}, \typ
317{rounded}, \typ {toposit} and \typ {fromposit}. The usual operators are also
318supported: \type{}, \type{}, \type{*}, \type{}, \type{}, as well as the
319binary \type {}. \type {}, \type {}, \type {<<} and \type {>>}. We can compare
320with \type {==}, \type {>=}, \type {<=} and \type {=}. The more verbose \type
321{bor}, \type {bxor}, \type {band}, \type {shift}, \type {rotate} are there too.
322
323There is a subset of math provided: \type {min}, \type {max}, \type {abs}, \type
324{conj}, \type {modf}, \type {acos}, \type {asin}, \type {atan}, \type {ceil},
325\type {cos}, \type {exp}, \type {exp2}, \type {floor}, \type {log}, \type
326{log10}, \type {log1p}, \type {log2}, \type {logb}, \type {pow}, \type {round},
327\type {sin}, \type {sqrt} and \type {tan}. Somewhat special are \type {NaN} and
328\type {NaR}.
329
330Currently integer division (\type {}) and modulo (\type {
331but that might happen at some time.
332
333\stopsectionlevel
334
335\stopdocument
336 |