-
Notifications
You must be signed in to change notification settings - Fork 0
/
14_trial_runs.html
744 lines (711 loc) · 55.5 KB
/
14_trial_runs.html
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
<!DOCTYPE html>
<html lang="en">
<head>
<meta charset="utf-8" />
<meta name="viewport" content="width=device-width, initial-scale=1.0" />
<meta charset="utf-8">
<meta name="viewport" content="width=device-width, initial-scale=1, shrink-to-fit=no">
<meta http-equiv="x-ua-compatible" content="ie=edge">
<title>14 试运行：Unums 面临计算挑战 — THE END of ERROR - Unum Computing 0.1 documentation</title>
<link rel="stylesheet" href="_static/material-design-lite-1.3.0/material.blue-deep_orange.min.css" type="text/css" />
<link rel="stylesheet" href="_static/sphinx_materialdesign_theme.css" type="text/css" />
<link rel="stylesheet" href="_static/fontawesome/all.css" type="text/css" />
<link rel="stylesheet" href="_static/fonts.css" type="text/css" />
<link rel="stylesheet" type="text/css" href="_static/pygments.css" />
<link rel="stylesheet" type="text/css" href="_static/basic.css" />
<link rel="stylesheet" type="text/css" href="_static/d2l.css" />
<script data-url_root="./" id="documentation_options" src="_static/documentation_options.js"></script>
<script src="_static/jquery.js"></script>
<script src="_static/underscore.js"></script>
<script src="_static/_sphinx_javascript_frameworks_compat.js"></script>
<script src="_static/doctools.js"></script>
<script src="_static/sphinx_highlight.js"></script>
<script src="_static/d2l.js"></script>
<script async="async" src="https://cdn.jsdelivr.net/npm/mathjax@3/es5/tex-mml-chtml.js"></script>
<link rel="index" title="Index" href="genindex.html" />
<link rel="search" title="Search" href="search.html" />
<link rel="next" title="小结" href="part1_summary.html" />
<link rel="prev" title="13 融合操作（一次性表达式）" href="13_fused_operations.html" />
</head>
<body>
<div class="mdl-layout mdl-js-layout mdl-layout--fixed-header mdl-layout--fixed-drawer"><header class="mdl-layout__header mdl-layout__header--waterfall ">
<div class="mdl-layout__header-row">
<nav class="mdl-navigation breadcrumb">
<a class="mdl-navigation__link is-active">14 试运行：Unums 面临计算挑战</a>
</nav>
<div class="mdl-layout-spacer"></div>
<nav class="mdl-navigation">
<form class="form-inline pull-sm-right" action="search.html" method="get">
<div class="mdl-textfield mdl-js-textfield mdl-textfield--expandable mdl-textfield--floating-label mdl-textfield--align-right">
<label id="quick-search-icon" class="mdl-button mdl-js-button mdl-button--icon" for="waterfall-exp">
<i class="material-icons">search</i>
</label>
<div class="mdl-textfield__expandable-holder">
<input class="mdl-textfield__input" type="text" name="q" id="waterfall-exp" placeholder="Search" />
<input type="hidden" name="check_keywords" value="yes" />
<input type="hidden" name="area" value="default" />
</div>
</div>
<div class="mdl-tooltip" data-mdl-for="quick-search-icon">
Quick search
</div>
</form>
<a id="button-show-source"
class="mdl-button mdl-js-button mdl-button--icon"
href="_sources/14_trial_runs.rst.txt" rel="">
<i class="material-icons">code</i>
</a>
<div class="mdl-tooltip" data-mdl-for="button-show-source">
Show Source
</div>
</nav>
</div>
<div class="mdl-layout__header-row header-links">
<div class="mdl-layout-spacer"></div>
<nav class="mdl-navigation">
<a class="mdl-navigation__link" href="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/jszheng/TheEndOfError">
<i class="fab fa-github"></i>
Github
</a>
</nav>
</div>
</header><header class="mdl-layout__drawer">
<!-- Title -->
<span class="mdl-layout-title">
<a class="title" href="index.html">
<span class="title-text">
THE END of ERROR - Unum Computing
</span>
</a>
</span>
<div class="globaltoc">
<span class="mdl-layout-title toc">Table Of Contents</span>
<nav class="mdl-navigation">
<ul class="current">
<li class="toctree-l1"><a class="reference internal" href="Preface.html">Preface</a></li>
<li class="toctree-l1"><a class="reference internal" href="00_how_to_read.html">如何读这本书</a></li>
<li class="toctree-l1"><a class="reference internal" href="Part1.html">Part 1 一种新的数字格式Unum</a></li>
<li class="toctree-l1"><a class="reference internal" href="01_Overview.html">1 概论</a></li>
<li class="toctree-l1"><a class="reference internal" href="02_BuildUpUnumFormat.html">2. 构造unum的格式</a></li>
<li class="toctree-l1"><a class="reference internal" href="03_TheOriginalSin.html">3. 计算机算术的原罪</a></li>
<li class="toctree-l1"><a class="reference internal" href="04_unum_format.html">4. 完整的unum格式定义</a></li>
<li class="toctree-l1"><a class="reference internal" href="05_hidden_scratchpads_3_layers.html">5. 隐藏的草稿本和三个层次</a></li>
<li class="toctree-l1"><a class="reference internal" href="06_info_per_bit.html">6 每个比特的信息</a></li>
<li class="toctree-l1"><a class="reference internal" href="07_fixed_size_unum_storage.html">7 定长的unum存储</a></li>
<li class="toctree-l1"><a class="reference internal" href="08_comparison_operations.html">8 比较操作</a></li>
<li class="toctree-l1"><a class="reference internal" href="09_add_sub_unbias_round.html">9 加减法和无偏差舍入的迷</a></li>
<li class="toctree-l1"><a class="reference internal" href="10_mul_div.html">10 乘法和除法</a></li>
<li class="toctree-l1"><a class="reference internal" href="11_power.html">11 求幂</a></li>
<li class="toctree-l1"><a class="reference internal" href="12_other_important_unary_ops.html">12 其他重要的一元运算</a></li>
<li class="toctree-l1"><a class="reference internal" href="13_fused_operations.html">13 融合操作（一次性表达式）</a></li>
<li class="toctree-l1 current"><a class="current reference internal" href="#">14 试运行：Unums 面临计算挑战</a></li>
<li class="toctree-l1"><a class="reference internal" href="part1_summary.html">小结</a></li>
<li class="toctree-l1"><a class="reference internal" href="Part2.html">Part 2 - 一种新的解决方法 Ubox</a></li>
<li class="toctree-l1"><a class="reference internal" href="15_TheOtherKindOfError.html">15. 另外一种误差</a></li>
<li class="toctree-l1"><a class="reference internal" href="16_avoid_interval_arith_pitfalls.html">16 避免区间算术陷阱</a></li>
<li class="toctree-l1"><a class="reference internal" href="17_meaning_of_solve_equ.html">17 “解”方程到底是什么意思？</a></li>
<li class="toctree-l1"><a class="reference internal" href="18_permission_to_guess.html">18 准许猜测</a></li>
<li class="toctree-l1"><a class="reference internal" href="19_pendulums_done_correctly.html">19 摆的正确计算</a></li>
<li class="toctree-l1"><a class="reference internal" href="20_two_body_problem.html">20 二体问题(以及多体问题)</a></li>
<li class="toctree-l1"><a class="reference internal" href="21_calculus_evil.html">21 微积分被认为是邪恶的：离散物理</a></li>
<li class="toctree-l1"><a class="reference internal" href="22_end_of_error.html">22 错误的终结</a></li>
<li class="toctree-l1"><a class="reference internal" href="Glossary.html">词汇表</a></li>
</ul>
</nav>
</div>
</header>
<main class="mdl-layout__content" tabIndex="0">
<script type="text/javascript" src="_static/sphinx_materialdesign_theme.js "></script>
<header class="mdl-layout__drawer">
<!-- Title -->
<span class="mdl-layout-title">
<a class="title" href="index.html">
<span class="title-text">
THE END of ERROR - Unum Computing
</span>
</a>
</span>
<div class="globaltoc">
<span class="mdl-layout-title toc">Table Of Contents</span>
<nav class="mdl-navigation">
<ul class="current">
<li class="toctree-l1"><a class="reference internal" href="Preface.html">Preface</a></li>
<li class="toctree-l1"><a class="reference internal" href="00_how_to_read.html">如何读这本书</a></li>
<li class="toctree-l1"><a class="reference internal" href="Part1.html">Part 1 一种新的数字格式Unum</a></li>
<li class="toctree-l1"><a class="reference internal" href="01_Overview.html">1 概论</a></li>
<li class="toctree-l1"><a class="reference internal" href="02_BuildUpUnumFormat.html">2. 构造unum的格式</a></li>
<li class="toctree-l1"><a class="reference internal" href="03_TheOriginalSin.html">3. 计算机算术的原罪</a></li>
<li class="toctree-l1"><a class="reference internal" href="04_unum_format.html">4. 完整的unum格式定义</a></li>
<li class="toctree-l1"><a class="reference internal" href="05_hidden_scratchpads_3_layers.html">5. 隐藏的草稿本和三个层次</a></li>
<li class="toctree-l1"><a class="reference internal" href="06_info_per_bit.html">6 每个比特的信息</a></li>
<li class="toctree-l1"><a class="reference internal" href="07_fixed_size_unum_storage.html">7 定长的unum存储</a></li>
<li class="toctree-l1"><a class="reference internal" href="08_comparison_operations.html">8 比较操作</a></li>
<li class="toctree-l1"><a class="reference internal" href="09_add_sub_unbias_round.html">9 加减法和无偏差舍入的迷</a></li>
<li class="toctree-l1"><a class="reference internal" href="10_mul_div.html">10 乘法和除法</a></li>
<li class="toctree-l1"><a class="reference internal" href="11_power.html">11 求幂</a></li>
<li class="toctree-l1"><a class="reference internal" href="12_other_important_unary_ops.html">12 其他重要的一元运算</a></li>
<li class="toctree-l1"><a class="reference internal" href="13_fused_operations.html">13 融合操作（一次性表达式）</a></li>
<li class="toctree-l1 current"><a class="current reference internal" href="#">14 试运行：Unums 面临计算挑战</a></li>
<li class="toctree-l1"><a class="reference internal" href="part1_summary.html">小结</a></li>
<li class="toctree-l1"><a class="reference internal" href="Part2.html">Part 2 - 一种新的解决方法 Ubox</a></li>
<li class="toctree-l1"><a class="reference internal" href="15_TheOtherKindOfError.html">15. 另外一种误差</a></li>
<li class="toctree-l1"><a class="reference internal" href="16_avoid_interval_arith_pitfalls.html">16 避免区间算术陷阱</a></li>
<li class="toctree-l1"><a class="reference internal" href="17_meaning_of_solve_equ.html">17 “解”方程到底是什么意思？</a></li>
<li class="toctree-l1"><a class="reference internal" href="18_permission_to_guess.html">18 准许猜测</a></li>
<li class="toctree-l1"><a class="reference internal" href="19_pendulums_done_correctly.html">19 摆的正确计算</a></li>
<li class="toctree-l1"><a class="reference internal" href="20_two_body_problem.html">20 二体问题(以及多体问题)</a></li>
<li class="toctree-l1"><a class="reference internal" href="21_calculus_evil.html">21 微积分被认为是邪恶的：离散物理</a></li>
<li class="toctree-l1"><a class="reference internal" href="22_end_of_error.html">22 错误的终结</a></li>
<li class="toctree-l1"><a class="reference internal" href="Glossary.html">词汇表</a></li>
</ul>
</nav>
</div>
</header>
<div class="document">
<div class="page-content" role="main">
<div class="section" id="unums">
<h1>14 试运行：Unums 面临计算挑战<a class="headerlink" href="#unums" title="Permalink to this heading">¶</a></h1>
<div class="figure align-default" id="id3">
<img alt="_images/image-20230627165205705.png" src="_images/image-20230627165205705.png" />
<p class="caption"><span class="caption-number">Fig. 207 </span><span class="caption-text">image-20230627165205705</span><a class="headerlink" href="#id3" title="Permalink to this image">¶</a></p>
</div>
<blockquote>
<div><center><p>伯克利教授 William Kahan（威廉.卡汗），IEEE 754 标准的主要架构师。</p>
</center></div></blockquote>
<div class="section" id="ii">
<h2>14.1 浮点数II：卡汉之怒<a class="headerlink" href="#ii" title="Permalink to this heading">¶</a></h2>
<blockquote>
<div><p>标题<em>Floating point II: The Wrath of Kahan</em>是模仿<a class="reference external" href="https://en.wikipedia.org/wiki/Star_Trek_II:_The_Wrath_of_Khan">Star Trek II:
The Wrath of
Khan</a></p>
</div></blockquote>
<p>伯克利大学教授威廉·卡汉（William
Kahan）擅长寻找浮点数给出严重错误答案的数学问题的例子……以及提出的替代方案甚至比浮点数更糟糕的例子。
unums能承受Kahan的愤怒吗？ 这是卡汉的例子之一：迭代</p>
<div class="math notranslate nohighlight" id="equation-14-trial-runs-0">
<span class="eqno">(52)<a class="headerlink" href="#equation-14-trial-runs-0" title="Permalink to this equation">¶</a></span>\[U_{i+2} = 111 - \frac{1130}{u_{i+1}} + \frac{3000}{u_i u_{i+1}}\]</div>
<p>从 <span class="math notranslate nohighlight">\(u_0=2、u_1=4\)</span> 开始。 使用浮点数时，迭代开始收敛到
6，但随后转向 100，即使 6 是有效的稳态答案。 迭代对于 6
以上的数字是稳定的，但对于 6 以下的数字不稳定。如果数值计算完美，则
<span class="math notranslate nohighlight">\(u_i\)</span> 值永远不会低于 6；
然而，舍入误差导致它们陷入不稳定区域，然后径直走向
100。下一页的循环使用双精度浮点数显示了这一点：</p>
<div class="figure align-default" id="id4">
<img alt="_images/image-20230627170312618.png" src="_images/image-20230627170312618.png" />
<p class="caption"><span class="caption-number">Fig. 208 </span><span class="caption-text">image-20230627170312618</span><a class="headerlink" href="#id4" title="Permalink to this image">¶</a></p>
</div>
<div class="highlight-default notranslate"><div class="highlight"><pre><span></span><span class="mf">18.5</span>
<span class="mf">9.37838</span>
<span class="mf">7.80115</span>
<span class="mf">7.15441</span>
<span class="mf">6.80678</span>
<span class="mf">6.59263</span>
<span class="mf">6.44947</span>
<span class="mf">6.34845</span>
<span class="mf">6.27444</span>
<span class="mf">6.21868</span>
<span class="mf">6.17563</span>
<span class="mf">6.13895</span>
<span class="mf">6.06036</span>
<span class="mf">5.1783</span>
<span class="o">-</span><span class="mf">11.6233</span>
<span class="mf">158.376</span>
<span class="mf">102.235</span>
<span class="mf">100.132</span>
</pre></div>
</div>
<p>作为浮点替代方案的一个示例，请尝试使用双精度区间算术的迭代循环，该算法每个值消耗
128 位而不是 64 位。它们最好物有所值：</p>
<div class="figure align-default" id="id5">
<img alt="_images/image-20230627170446017.png" src="_images/image-20230627170446017.png" />
<p class="caption"><span class="caption-number">Fig. 209 </span><span class="caption-text">image-20230627170446017</span><a class="headerlink" href="#id5" title="Permalink to this image">¶</a></p>
</div>
<div class="figure align-default" id="id6">
<img alt="_images/image-20230627170532750.png" src="_images/image-20230627170532750.png" />
<p class="caption"><span class="caption-number">Fig. 210 </span><span class="caption-text">image-20230627170532750</span><a class="headerlink" href="#id6" title="Permalink to this image">¶</a></p>
</div>
<p>他们<strong>不</strong>值得。
这个答案是典型的区间算术：一个正确但无用的结果，表明答案是：<strong>实数线上的某个位置</strong>。
现在您可以开始明白为什么区间算术吸引了非常有限的粉丝俱乐部。</p>
<p>现在尝试使用
unums，调高精度（但不是指数，因为这些数字不是特别大或特别小）。
请原谅分数的长小数。 循环后的reportdatamotion 调用报告数据移动成本。</p>
<div class="figure align-default" id="id7">
<img alt="_images/image-20230627181152756.png" src="_images/image-20230627181152756.png" />
<p class="caption"><span class="caption-number">Fig. 211 </span><span class="caption-text">image-20230627181152756</span><a class="headerlink" href="#id7" title="Permalink to this image">¶</a></p>
</div>
<p>因此 unums 前往正确的数字 6。如果继续前进，左端点将降至 6
以下，然后向<span class="math notranslate nohighlight">\(-\infty\)</span>发散，但请记住，unum 环境可以使用
<strong>needmorefracQ</strong> 自动检测相对宽度何时变得过高，并停止 (
正如我们在这里所做的那样）或提高精度。 Unum
不会像浮点数那样表现出虚假收敛。 此外，它们使用比双精度区间算术（每个值
128 位）更少的位来完成所有这些操作。</p>
<p>许多科学程序员普遍认为以下神话是防范浮点错误危险的简单方法：</p>
<table border="2"><tr><td bgcolor="lightyellow"><p>如果您不断增加浮点数的精度，并且答案是一致的，那么答案对于该精度来说是正确的。</p>
</td></tr></table><p>卡汉构建了以下极其聪明的例子来证明这个神话是多么错误。
对于浮点数，即使精度非常高，您也经常得到 0；</p>
<p>正确答案其实是 1。</p>
<p>定义函数：</p>
<ul class="simple">
<li><p><span class="math notranslate nohighlight">\(E(0)=1\)</span></p></li>
<li><p><span class="math notranslate nohighlight">\(E(z)=\frac{e^z-1}{z}\)</span></p></li>
<li><p><span class="math notranslate nohighlight">\(Q[x]=\left|x-\sqrt{x^2+1}\right|-\frac{1}{x+\sqrt{x^2+1}}\)</span></p></li>
<li><p><span class="math notranslate nohighlight">\(H(x)=E(Q(X))^2\)</span></p></li>
</ul>
<p>计算<span class="math notranslate nohighlight">\(x=15.0, 16.0, 17.0, 9999.0\)</span>下的<span class="math notranslate nohighlight">\(H(x)\)</span>的值，换个精度，比如用BigDecimal再做一次看看</p>
<p>“BigDecimal”是 Java
中的一种数据类型，允许扩展精度算术，通过使用大量存储，在用户可见层中允许更多数字来处理舍入误差，这是历史上的众多尝试之一。
卡汉教授确切地知道如何揭露浮点数学的缺陷。 使用双精度浮点数，以下是您对
x 的所有四个建议值得到的错误答案：</p>
<div class="figure align-default" id="id8">
<img alt="_images/image-20230627195316582.png" src="_images/image-20230627195316582.png" />
<p class="caption"><span class="caption-number">Fig. 212 </span><span class="caption-text">image-20230627195316582</span><a class="headerlink" href="#id8" title="Permalink to this image">¶</a></p>
</div>
<p>即使精度更高，例如 IEEE 四精度或 BigDecimal
扩展精度，它几乎每次仍然会出现错误。 x 的一些分散值给出了正确答案：{1,
1, 1, 1}，可以通过符号而不是数字来计算代数。
关键点是，当答案完全错误时，简单地提高精度并得到相同、完全一致的答案就可能发生。</p>
<p>那么传统的区间运算呢？它应该是“严格的”？
为了给传统区间提供所有可能的优势，我们甚至可以让 Mathematica
以符号方式计算结果，并仅在最后使用数值浮点数对其进行求值，方法是重新定义
e、q 和 h 函数以使用区间算术。 区间为零的测试与我们使用 unums
测试它的方式类似，其中如果零在区间 z 中，则
<strong>IntervalMemberQ</strong>[z,<strong>Interval</strong>[{0,0}]] 为 True，否则为 False：</p>
<div class="figure align-default" id="id9">
<img alt="_images/image-20230627195919811.png" src="_images/image-20230627195919811.png" />
<p class="caption"><span class="caption-number">Fig. 213 </span><span class="caption-text">image-20230627195919811</span><a class="headerlink" href="#id9" title="Permalink to this image">¶</a></p>
</div>
<p>又一个史诗般的失败。
传统的区间算术再次为您提供整个实数轴作为可能的答案，然后迫使您分析出了什么问题并调试计算丢失所有信息的位置。
此外，与数值计算相比，符号数学的使用极其昂贵且缓慢。
这个例子看起来是一个很难破解的难题！</p>
<table border="2"><tr><td bgcolor="lightblue"><p>像本书这样的“戏剧性”时刻并不多，但如果有的话，那就是它了。 unums
能打败卡汉教授的怪物吗？ 也许通过使用一些非常高的环境设置来提高精度。</p>
</td></tr></table><p>这是函数定义的 unum 版本：</p>
<div class="figure align-default" id="id10">
<img alt="_images/image-20230627200307793.png" src="_images/image-20230627200307793.png" />
<p class="caption"><span class="caption-number">Fig. 214 </span><span class="caption-text">image-20230627200307793</span><a class="headerlink" href="#id10" title="Permalink to this image">¶</a></p>
</div>
<p>只是为了好玩，我们展示了其中最小、最朴素的 unum，即四位 Warlpiri
number集。 在这场战斗中抱有成功希望的唯一原因是，所有
unum，即使是那些只有一位分数和一位指数的
unum，都正在包装一个强大的武器，这是其他计算数字系统所不具备的：绝对诚实
他们对某个值做了什么和不知道什么。下面就是了</p>
<div class="figure align-default" id="id11">
<img alt="_images/image-20230627200445542.png" src="_images/image-20230627200445542.png" />
<p class="caption"><span class="caption-number">Fig. 215 </span><span class="caption-text">image-20230627200445542</span><a class="headerlink" href="#id11" title="Permalink to this image">¶</a></p>
</div>
<p>Warlpiri 的数学是正确的，尽管它不能数到超过二！ Unums
用每个数字不到七位的存储空间计算出正确的答案。</p>
<div class="figure align-default" id="id12">
<img alt="_images/image-20230627200608612.png" src="_images/image-20230627200608612.png" />
<p class="caption"><span class="caption-number">Fig. 216 </span><span class="caption-text">image-20230627200608612</span><a class="headerlink" href="#id12" title="Permalink to this image">¶</a></p>
</div>
</div>
<div class="section" id="rump">
<h2>14.2 Rump的烦恼事<a class="headerlink" href="#rump" title="Permalink to this heading">¶</a></h2>
<p>卡汉并不是唯一一个想出办法揭露计算机数字系统缺点的人。 IBM 科学家
Siegfried Rump
提出了一个很好的例子，说明为什么仅仅尝试不同的浮点精度并寻求一致性是不够的</p>
<p>求：</p>
<p><span class="math notranslate nohighlight">\(333.75 y^6 + x^2(11 x^2 y^2 - y^6 -121 y^4 -2)+5.5y^8 + \frac{x}{2y}\)</span></p>
<p>其中 x=77 617, y=33 096</p>
<p>以下是 Rump 的公式和输入参数，以及现代 IEEE 双精度浮点数告诉我们的内容</p>
<div class="figure align-default" id="id13">
<img alt="_images/image-20230627201112787.png" src="_images/image-20230627201112787.png" />
<p class="caption"><span class="caption-number">Fig. 217 </span><span class="caption-text">image-20230627201112787</span><a class="headerlink" href="#id13" title="Permalink to this image">¶</a></p>
</div>
<p>对于 IEEE 风格的浮点数，Rump
的公式对于单精度、双精度甚至四精度给出了相同的答案，约为
<span class="math notranslate nohighlight">\(1.18\times10^{21}\)</span>，但在这三种情况下都存在严重错误。 正确答案是
<span class="math notranslate nohighlight">\(–0.827\cdots\)</span>，因此浮点数甚至无法得到正确答案的符号。
更不用说相差 21 个数量级的小问题了。</p>
<p>Rump 在 1970 年代老式 IBM 大型机上写程序得到该函数的计算结果为 1.172603
单精度， 1.1726039400531双精度 ，以及 1.172603940053178
IBM“扩展”（我们现在称为四元）精度中的 。</p>
<p>IBM 370 浮点格式使用以 16 为基数的指数，并且与现代 IEEE
浮点有其他差异，这就是为什么它会产生如此不同的三个高度一致的结果。
任何试图检查计算可靠性的人都会得出结论：结果是有效的，并且只需要单精度。
请注意，公式中使用的每个数字都可以表示为精确浮点数：333.75, 11, 121, 2,
5.5，并且两个整数输入都是除以 2 的幂的整数，因此在此
情况下，简单地输入问题不会引入任何舍入误差。</p>
<p>如果不是以精确结果的形式呈现，那么偏差 21
个数量级以及符号错误可能会更容易被原谅。
区间算术可以来拯救并保证包含正确的结果吗？ 或许吧。</p>
<p>尝试使用名为 rumpint 的rump版本，它使用 128 位间隔（即每个端点使用 64
位浮点数），至少您会得到答案的界限：</p>
<div class="figure align-default" id="id14">
<img alt="_images/image-20230627202131622.png" src="_images/image-20230627202131622.png" />
<p class="caption"><span class="caption-number">Fig. 218 </span><span class="caption-text">image-20230627202131622</span><a class="headerlink" href="#id14" title="Permalink to this image">¶</a></p>
</div>
<p>一个巨大的间隔，但也是一个界限。
至少传统的间隔表明rump问题可能容易出现严重错误。
（不幸的是，传统区间也经常暗示着完全稳定的计算。）像往常一样，传统区间结果的巨大宽度让我们对下一步该做什么感到困惑。
尝试使用四精度边界的传统间隔？ 他们也失败了。</p>
<p>用 unum 试试。 因为关于 unums 正在测试的一件事是易用性，所以这里是对
unum 版本的以下代码的一些解释。
<span class="math notranslate nohighlight">\(\otimes\)</span>运算符优先于<span class="math notranslate nohighlight">\(\oplus\)</span>和<span class="math notranslate nohighlight">\(\ominus\)</span>
运算符，就像传统公式中乘法优先于加法和减法一样，因此我们不必使用括号来编写，例如
<span class="math notranslate nohighlight">\(a \otimes b \oplus c \otimes d\)</span>。 然而，为了“打破
<span class="math notranslate nohighlight">\(x \oplus y \oplus z\)</span>之类的表达式的束缚”，我们使用括号从左到右对其进行分组：<span class="math notranslate nohighlight">\((x \oplus y) \oplus z\)</span>。
在 <span class="math notranslate nohighlight">\(\oplus\)</span> 和 <span class="math notranslate nohighlight">\(\ominus\)</span>
运算符周围留出足够的空间也有助于提高可读性。 为了清楚起见，我们为
<span class="math notranslate nohighlight">\(x^2、y^2、y^4、y^6 和 y^8\)</span> 引入一些临时变量。 <strong>rumpu</strong> 函数使用
ubounds 计算 Rump 的公式：</p>
<div class="figure align-default" id="id15">
<img alt="_images/image-20230627202926314.png" src="_images/image-20230627202926314.png" />
<p class="caption"><span class="caption-number">Fig. 219 </span><span class="caption-text">image-20230627202926314</span><a class="headerlink" href="#id15" title="Permalink to this image">¶</a></p>
</div>
<p>假设我们希望答案的精度至少达到两位小数，因此我们设置 relwidthtolerance =
0.005。 从 Warlpiri 环境 {0, 0}
开始，让计算机自动查找何时有足够的位来获得可接受的答案，如第 9.3
节所示的流程图。 事实证明，{3, 7}
是获得可接受精度的环境，并且它确实比两位小数要好得多：</p>
<div class="figure align-default" id="id16">
<img alt="_images/image-20230627203053385.png" src="_images/image-20230627203053385.png" />
<p class="caption"><span class="caption-number">Fig. 220 </span><span class="caption-text">image-20230627203053385</span><a class="headerlink" href="#id16" title="Permalink to this image">¶</a></p>
</div>
<p>unum 方法将答案的精确度限制在 39 位小数，并且每个数字平均使用 70 位。
unums 在这里表现出色的主要原因是它们可以根据需要动态调整精度。
创建表达式的中间结果非常大，需要 37 位小数的精度才能准确表示。
这些小数稍后会在计算中被抵消，从而产生 unum 格式用很少位存储的小整数。
携带 11 位 utag 的开销会带来数倍的回报。
如果我们使用融合运算，每个数字的位数会更低。</p>
<table border="2"><tr><td bgcolor="lightblue"><center><p>更少的比特，更好的答案</p>
</center></td></tr></table></div>
<div class="section" id="id1">
<h2>14.3 二次公式<a class="headerlink" href="#id1" title="Permalink to this heading">¶</a></h2>
<p>前面的示例经过精心设计，旨在暴露浮动问题。
读者可能会认为这些示例过于做作，无法反映大多数计算任务，而浮点数似乎做得相当好。
因此，下一个示例展示了一个来自基本主流计算的示例，其中舍入错误会给粗心的程序员带来不愉快的意外。
恰好有一个处理二次公式的特定问题的技巧，并且数值分析书籍热切地告诉您如何调用该技巧，但我们并不是在寻找程序员必须学习的技巧的集合;
计算机应该尽可能多地自动化数值分析。
即使以简单的方式，也就是浮点数通常使用的方式，使用 Unum 也应该可以工作。</p>
<p>二次方程 <span class="math notranslate nohighlight">\(a x^2 + b x + c = 0\)</span> 的根 <span class="math notranslate nohighlight">\(r_1\)</span> 和 <span class="math notranslate nohighlight">\(r_2\)</span>
的公式具有数值风险，我们可以使用它来比较 float 与
unum，无论是为了准确性还是为了存储效率：</p>
<div class="math notranslate nohighlight" id="equation-14-trial-runs-1">
<span class="eqno">(53)<a class="headerlink" href="#equation-14-trial-runs-1" title="Permalink to this equation">¶</a></span>\[r_{1,2}=\frac{-b\pm\sqrt{b^2-4ac}}{2a}\]</div>
<p>当然，平方根通常会产生不精确的结果。 更微妙的危险是，当 <span class="math notranslate nohighlight">\(b^2\)</span>
相对于 <span class="math notranslate nohighlight">\(4 a c\)</span> 太大时，<span class="math notranslate nohighlight">\(\sqrt{b^2 - 4 a c}\)</span> 将与 <span class="math notranslate nohighlight">\(b\)</span>
共享许多前导数字。 这会导致计算 <span class="math notranslate nohighlight">\(-b \pm \sqrt{b2 - 4 ac}\)</span>
时出现数字抵消，因为 ± 的 + 或 –
情况都会导致类似数字相减，只留下几个有效数字。</p>
<p>对于 a = 3、b = 100 和 c = 2，请注意有一个根与零有多接近：</p>
<div class="figure align-default" id="id17">
<img alt="_images/image-20230627205308483.png" src="_images/image-20230627205308483.png" />
<p class="caption"><span class="caption-number">Fig. 221 </span><span class="caption-text">image-20230627205308483</span><a class="headerlink" href="#id17" title="Permalink to this image">¶</a></p>
</div>
<p>以下是两个根的精确表达式，其值精确到小数点后 7 位，供参考：</p>
<div class="figure align-default" id="id18">
<img alt="_images/image-20230627205515718.png" src="_images/image-20230627205515718.png" />
<p class="caption"><span class="caption-number">Fig. 222 </span><span class="caption-text">image-20230627205515718</span><a class="headerlink" href="#id18" title="Permalink to this image">¶</a></p>
</div>
<p>在IEEE单精度中，b2=100可以精确表示； 4 a c = 24 和 b2 - 4 a c = 9976
也是如此。到目前为止，一切都很好。 9976 的平方根是： 99.8799279…</p>
<p>但是单精度浮点数要求我们接受以下不正确的数字，我们在正确的数字之外使用橙色数字，并在末尾使用“<span class="math notranslate nohighlight">\(\downarrow\)</span>”符号来强调显示的小数扩展正是该浮点数所表达的精确值：</p>
<div class="figure align-default" id="id19">
<img alt="_images/image-20230627205811300.png" src="_images/image-20230627205811300.png" />
<p class="caption"><span class="caption-number">Fig. 223 </span><span class="caption-text">image-20230627205811300</span><a class="headerlink" href="#id19" title="Permalink to this image">¶</a></p>
</div>
<p>这肯定在正确平方根的 1 个 ULP 范围内。 问题出现在下一步的计算中，将-b =
-100 添加到上面的数字中。
这样做会消除左边的两个数字，并且错误的橙色数字会放大相对误差：</p>
<div class="figure align-default" id="id20">
<img alt="_images/image-20230627210153376.png" src="_images/image-20230627210153376.png" />
<p class="caption"><span class="caption-number">Fig. 224 </span><span class="caption-text">image-20230627210153376</span><a class="headerlink" href="#id20" title="Permalink to this image">¶</a></p>
</div>
<p>在最后一步中，除以 2a，我们最终与参考答案只有四位小数一致，结果错误超过
100 个 ULP：</p>
<div class="figure align-default" id="id21">
<img alt="_images/image-20230627210353400.png" src="_images/image-20230627210353400.png" />
<p class="caption"><span class="caption-number">Fig. 225 </span><span class="caption-text">image-20230627210353400</span><a class="headerlink" href="#id21" title="Permalink to this image">¶</a></p>
</div>
<p>如果计算机跟踪了 ULP
误差是如何放大的，它就会知道将答案表示为“~-0.02001”或只有四位小数的数字。
相反，作为单精度答案，它被欺骗性地显示为“–0.02001190”，没有任何信息让我们知道最后三个看起来精确的小数是错误的。</p>
<p>数值分析书籍中教授的技巧是使用恒等式 r1 r2 = c/a 来求 r1。 由于 r2
不会经历数字抵消，因此精度通常保持不变。
然而，大多数人在学习二次公式后很快就会忘记它的特殊性，并且无论如何，要求程序员记住并应用这样的技巧似乎是不必要的负担。
相反，尝试使用 unum 环境，这是可以模拟单精度浮点数的最小环境：</p>
<div class="highlight-mathematica notranslate"><div class="highlight"><pre><span></span><span class="n">setenv</span><span class="p">[{</span><span class="mi">3</span><span class="p">,</span><span class="w"> </span><span class="mi">5</span><span class="p">}]</span><span class="w"></span>
</pre></div>
</div>
<p>这是接近零的根的二次公式的 unum 版本，使用了之前开发的算术运算。
重置移动位数和个数计数器，以便我们可以将每个值的位数与单精度所需的 32
位进行比较：</p>
<div class="figure align-default" id="id22">
<img alt="_images/image-20230627211047288.png" src="_images/image-20230627211047288.png" />
<p class="caption"><span class="caption-number">Fig. 226 </span><span class="caption-text">image-20230627211047288</span><a class="headerlink" href="#id22" title="Permalink to this image">¶</a></p>
</div>
<p>界限精确到小数点后七位：–0.02001201… {3, 5} unum
环境计算答案中的七个有效小数并严格限制错误，但每个数字使用的平均位数少于
25 位。 因此，在这种情况下，尽管 utag 存在开销，但 unum
提供的精度几乎是标准 IEEE 浮点数的两倍，但使用的位数更少。
当有数字cancel现象时，Unum 存储空间会自动缩小。
如果我们有使用<strong>unify</strong>操作的策略，我们可以自动使这变得更加经济。
最终的 87 位 ubound 答案统一为 41 位，即一个 unum：</p>
<div class="figure align-default" id="id23">
<img alt="_images/image-20230627211428316.png" src="_images/image-20230627211428316.png" />
<p class="caption"><span class="caption-number">Fig. 227 </span><span class="caption-text">image-20230627211428316</span><a class="headerlink" href="#id23" title="Permalink to this image">¶</a></p>
</div>
<div class="figure align-default" id="id24">
<img alt="_images/image-20230627211507319.png" src="_images/image-20230627211507319.png" />
<p class="caption"><span class="caption-number">Fig. 228 </span><span class="caption-text">image-20230627211507319</span><a class="headerlink" href="#id24" title="Permalink to this image">¶</a></p>
</div>
<p><strong>unify</strong> 操作仅将边界宽度扩展了约 1.5 倍，同时减少了约 2.1
倍的存储空间，因此在这种情况下，<strong>unify</strong> 增加了每比特的信息量。
结果仍然好到小数点后七位。</p>
<p>虽然我们没有在这里展示，但如果将 {3, 6} 环境结果与双精度 IEEE
结果进行比较，unum 会得到 16 位有效数字，而使用双精度则为 12 位，平均
unum 大小为 33 位，而双精度则为 64 位 。
我们<em>以大约一半的位数获得了超过 33% 的准确度</em>。</p>
<p>这个例子的要点是：<strong>不需要任何技巧</strong>。
二次公式示例支持这样的主张：除了比浮点数更准确且使用更少的存储空间之外，unum
还比浮点数更易于使用。 与浮点数相比，unum
的计算更有可能适应程序员想要计算的任何公式，而无需使用代数恒等式对计算进行复杂的重新排列。</p>
<table border="2"><tr><td bgcolor="lightblue"><center><p>“我们有如此多的数学技巧，以至于我们永远不会停止告诉人们如何做事……物理老师总是展示技术，而不是如何解决物理问题的精神。”
——理查德·费曼</p>
</center></td></tr></table></div>
<div class="section" id="id2">
<h2>14.4 贝利的数字噩梦<a class="headerlink" href="#id2" title="Permalink to this heading">¶</a></h2>
<p>David Bailey
是一位数值计算专家，他建立了扩展精度算术库，他提出了一个非常不稳定的系统，该系统由两个未知数的两个方程组成，用作测试：</p>
<div class="math notranslate nohighlight" id="equation-14-trial-runs-2">
<span class="eqno">(54)<a class="headerlink" href="#equation-14-trial-runs-2" title="Permalink to this equation">¶</a></span>\[\begin{split}0.25510582 x + 0.52746197 y = 0.79981812 \\
0.80143857 x + 1.65707065 y = 2.51270273\end{split}\]</div>
<p>这些方程看起来实在是太简单了。
假设可以精确的输入系数的小数，该系统可以精确求解 x = –1 和 y = 2
。我们可以尝试高斯消去法或克莱默法则。
对于这么小的系统，克莱默规则实际上看起来更简单。 回想一下克莱默规则：</p>
<div class="math notranslate nohighlight" id="equation-14-trial-runs-3">
<span class="eqno">(55)<a class="headerlink" href="#equation-14-trial-runs-3" title="Permalink to this equation">¶</a></span>\[\begin{split}\text{解方程 } \left ( \begin{array}{l} ax+by=u \\ cx + dy = v \end{array} \right )
\\ \text{对于x和y，计算分母。 } D=ad-bc \\ \text{ 如果D不为0，那么 } x=\frac{ud-vb}{D} \text{ 且 } y=\frac{av-cu}{D}\end{split}\]</div>
<p>从几何角度来说，求解两个未知数的两个方程意味着找到两条线的交点：</p>
<div class="figure align-default" id="id25">
<img alt="_images/image-20230628082834650.png" src="_images/image-20230628082834650.png" />
<p class="caption"><span class="caption-number">Fig. 229 </span><span class="caption-text">image-20230628082834650</span><a class="headerlink" href="#id25" title="Permalink to this image">¶</a></p>
</div>
<p>行列式D越小，直线越接近平行。 如果 D
恰好为零，则这些线要么是分开的并且从不相交，要么它们可能是同一条线并且到处相交。
假设线条几乎平行，并且输入值 a、b、c、d、u 和 v 中存在一点“摆动”，这是由
±0.5 ULP 的舍入误差引起的。
交点就会到处移动，这意味着问题不稳定：输入的微小变化会导致输出的巨大变化。
在贝利噩梦般的不稳定问题中，精确的确定 D 是</p>
<div class="math notranslate nohighlight" id="equation-14-trial-runs-4">
<span class="eqno">(56)<a class="headerlink" href="#equation-14-trial-runs-4" title="Permalink to this equation">¶</a></span>\[\begin{split}D = 0.25510582 \times 1.65707065 - 0.52746197 \times 0.80143857 \\
= 0.4227283669661830 - 0.4227283669661829 \\
= 0.0000000000000001 = 10^{-16}\end{split}\]</div>
<p>行列式不为零，但它很小，即使输入中 0.5 ULP 摆动也会导致答案与答案 x =
-1, y = 2 发生约 ±1 亿的变化。行列式很小没有问题 只要是准确的；
但因为这里是通过减去非常相似的浮点值来计算的，因此计算中的不准确性被放大，包括从十进制输入到二进制格式的转换。</p>
<p>由于 IEEE 双精度能够达到几乎 16 位小数精度，而输入只有 8 或 9
位小数，因此它似乎可以做到这一点。 <strong>但是它失败了</strong>。
当您使用双精度和克莱默规则时，减法会导致大量有效数字丢失。
行列式应该正好是 10-16：</p>
<div class="math notranslate nohighlight" id="equation-14-trial-runs-5">
<span class="eqno">(57)<a class="headerlink" href="#equation-14-trial-runs-5" title="Permalink to this equation">¶</a></span>\[det = 0.25510582 \times 1.65707065 - 0.52746197 \times 0.80143857\]</div>
<p><span class="math notranslate nohighlight">\(1.66533 \times 10^{-16}\)</span></p>
<p>行列式不为零，因此通过克莱默法则计算 x 和 y：</p>
<div class="math notranslate nohighlight" id="equation-14-trial-runs-6">
<span class="eqno">(58)<a class="headerlink" href="#equation-14-trial-runs-6" title="Permalink to this equation">¶</a></span>\[\begin{split}{\begin{array}{l}
(0.79981812 \times 1.65707065 - 2.51270273 \times .52746197) \div det \\
(0.25510582 \times 2.51270273 - .80143857 \times 0.79981812) \div det
\end{array}}\end{split}\]</div>
<p><span class="math notranslate nohighlight">\(\{ 0., 1.33333 \}\)</span></p>
<p>该答案与正确答案 {x, y} = {-1, 2} 没有任何相似之处。
至此，读者不会感到非常惊讶，尽管使用两倍的位数来表示数字，但双精度区间算术也会失败。
x 和 y 的计算结果均为 <span class="math notranslate nohighlight">\([-\infty, \infty]\)</span>（采用传统间隔）</p>
<p>与所有输入值都可以用浮点数精确表达的 Rump 问题不同，Bailey
问题有两个错误来源：将十进制输入转换为二进制浮点数，以及中间结果的精度不足。
我们可以将这两种效果分开。
要求程序员知道这个理论知识并不过分：那就是任何线性系统可以进行缩放而答案不改变。
由于浮点数可以精确地表示整数范围，因此只需将问题缩放 <span class="math notranslate nohighlight">\(10^8\)</span>
以使所有值都为整数，然后重试：</p>
<div class="math notranslate nohighlight" id="equation-14-trial-runs-7">
<span class="eqno">(59)<a class="headerlink" href="#equation-14-trial-runs-7" title="Permalink to this equation">¶</a></span>\[det = 25 510 582. \times 165 707 065. - 52 746 197. \times 80 143 857.\]</div>
<p><span class="math notranslate nohighlight">\(1.\)</span></p>
<p>这次行列式被精确地计算出来了。 这样就消除了小数到浮点转换的问题。
所以，砥砺前行：</p>
<div class="figure align-default" id="id26">
<img alt="_images/image-20230628085625387.png" src="_images/image-20230628085625387.png" />
<p class="caption"><span class="caption-number">Fig. 230 </span><span class="caption-text">image-20230628085625387</span><a class="headerlink" href="#id26" title="Permalink to this image">¶</a></p>
</div>
<p>x 仍然是错的，但至少 y 是对的。 x
中的误差是由于第二个缺点，即中间结果的精度不够。 问题是浮点仅表达精确到
ULP 变为宽度 2 而不是宽度 1 的大小的整数。第一个乘法
<span class="math notranslate nohighlight">\(79981812 \times 165707065\)</span> 应该是</p>
<blockquote>
<div><center><p>13 253 551 319 901 780</p>
</center></div></blockquote>
<p>这恰好可以用双精度表示，因为它是偶数。而乘积<span class="math notranslate nohighlight">\(251270273 \times 52746197\)</span>结果应该是</p>
<blockquote>
<div><center><p>13 253 551 319 901 781</p>
</center></div></blockquote>
<p>但该整数无法表达，因此它会向下舍入到最接近的偶数，从而使两个乘积之间的差值为
0 而不是 -1。</p>
<p>Bailey 认为供应商应该开始在硬件中支持四精度 IEEE
浮点，以便更容易解决此类问题。 Unums
指出了一种不同的方法，这种方法不要求程序员为了得到像 –1 和 2
这样的答案而使用 34 位十进制值：而是要<strong>依靠暂存器</strong>。 像
<span class="math notranslate nohighlight">\(u d - v b\)</span>
这样的克莱默规则计算的值只是点积，并且融合点积仅在所有乘积累加后才进行舍入。
不要采用四精度，而是采用另一种方式：使用 {3, 5}, 环境，可以模拟 IEEE
单精度。 输入值是可以精确表示的。
通过在原型中使用融合点积，程序员不必担心中间计算占用更多空间。</p>
<div class="figure align-default" id="id27">
<img alt="_images/image-20230628090316194.png" src="_images/image-20230628090316194.png" />
<p class="caption"><span class="caption-number">Fig. 231 </span><span class="caption-text">image-20230628090316194</span><a class="headerlink" href="#id27" title="Permalink to this image">¶</a></p>
</div>
<p>这就可以了，用了不到不准确的双精度解决方案所使用的一半的位数，就得到了计算准确的答案。
讽刺的是，现在可以使用最原始的 unum 环境来存储 Bailey
的苛刻问题的答案，因为全零 utag 左侧的 101 始终代表 –1，而全零 utag
左侧的 010 始终代表 2。</p>
<p>这里的教训是，所选环境至少最初应与输入值精度和要求输出值的精度相匹配，而
g-layer可以处理临时更精确的临时计算的需要。 在第 2
部分中，我们将探索一种技术，当输入值是范围而不是精确值时，该技术可以自动找到两个未知数中的两个方程问题的完整解集。
当问题不适定时，就像贝利的问题一样，这种技术会可以立即将其暴露出来。</p>
</div>
<div class="section" id="unum">
<h2>14.5 使用unum做快速傅里叶变换<a class="headerlink" href="#unum" title="Permalink to this heading">¶</a></h2>
<div class="figure align-default" id="id28">
<img alt="_images/image-20230628092444476.png" src="_images/image-20230628092444476.png" />
<p class="caption"><span class="caption-number">Fig. 232 </span><span class="caption-text">image-20230628092444476</span><a class="headerlink" href="#id28" title="Permalink to this image">¶</a></p>
</div>
<blockquote>
<div><p>石油和天然气矿床的勘探是通过巨大的地震探测器阵列来完成的，这些探测器记录强声脉冲的回声。
这可能是最古老的“大数据”问题，因为记录下来进行处理的比特还真是通过卡车运回的。
勘探行业长期以来一直在努力使数据移动和存储更加经济。 Unum
是做到这一点的关键，特别是可以减少最重要的内核信号处理的数据移动成本：快速傅里叶变换。</p>
</div></blockquote>
<p>本节涉及一些相当高阶的数学，因此请将其视为可选读。
信号处理是计算机最重要的用途之一，也是 unum 可以发挥特别作用的领域。
信号处理硬件在消费电子产品中无处不在，它可以改善相机和电视图像，为音频播放器编码和播放声音，甚至控制汽车系统。
数字信号处理以及一般科学和工程中最重要的运算之一是傅里叶变换。
如果您听到一个音乐和弦并且有足够好的耳朵来辨别各个音符，那么您所做的就是傅立叶变换的作用。
它将信号或任何函数转换为不同频率和幅度的正弦和余弦函数。
傅里叶逆变换将频率分量表示转换回原始信号。</p>
<p>计算离散输入数据的傅里叶变换的常用方法是所谓的快速傅里叶变换（简称
FFT）的某种变体。 对于具有 n 个数据点的信号，FFT 算法会通过大约 $5 n
log_2(n) $次浮点运算生成傅立叶变换。
在现代计算机上，性能受到严重限制的不是算术量，而是移动数据的需要。
如果我们可以安全地减少操作数中的位数，例如使用
unum，我们可以按比例运行得更快，而且我们还可以得到正确的答案范围。</p>
<p>许多信号信息来自仅精确到几个位的模数转换器。 目前声音编码的标准是每通道
16 位精度，但对于许多应用（如石油勘探）来说，数据只有 12
位精度（相当于大约 3.6 位十进制数字）。
如果您查看音响设备上的图形均衡器（或音乐软件显示的图形均衡器），您会看到傅里叶变换的作用。
它们连续实时计算从最低（低音）到最高（高音）的频率强度：</p>
<div class="figure align-default" id="id29">
<img alt="_images/image-20230628093502627.png" src="_images/image-20230628093502627.png" />
<p class="caption"><span class="caption-number">Fig. 233 </span><span class="caption-text">image-20230628093502627</span><a class="headerlink" href="#id29" title="Permalink to this image">¶</a></p>
</div>
<p>假设信号范围从 –8 到 8，并且只有 12 或 16 位小数。
在这种情况下，我们需要的只是一个小的 unum
标签来表示模数转换器可以提供的所有内容：</p>
<div class="figure align-default" id="id30">
<img alt="_images/image-20230628093543465.png" src="_images/image-20230628093543465.png" />
<p class="caption"><span class="caption-number">Fig. 234 </span><span class="caption-text">image-20230628093543465</span><a class="headerlink" href="#id30" title="Permalink to this image">¶</a></p>
</div>
<p>信号处理的标准基准测试长期以来一直是“1K CFFT”，即 1024
点复快速傅里叶变换。 下一页是该算法的简单通用长度版本，用于比较浮点数与
unum。 像 gg 和 ww 这样的双字母变量是 <span class="math notranslate nohighlight">\(a + bi\)</span> 形式的复数。</p>
<div class="figure align-default" id="id31">
<img alt="_images/image-20230628093724725.png" src="_images/image-20230628093724725.png" />
<p class="caption"><span class="caption-number">Fig. 235 </span><span class="caption-text">image-20230628093724725</span><a class="headerlink" href="#id31" title="Permalink to this image">¶</a></p>
</div>
<p>如果输入数据只有 12 到 16 个有效位，为什么不简单地使用 IEEE
半精度浮点数来节省存储空间呢？ 原因之一是 IEEE 半精度有 10
位小数，这还不够。 另一个原因是三角函数，对于 n 点 FFT，在 n
个位置对单位圆进行采样。 假设 n = 1024； <span class="math notranslate nohighlight">\(\frac{2 \pi}{1024}\)</span>
的余弦为 0.999981175 ，但是用10 位小数则四舍五入为数字 1。
<span class="math notranslate nohighlight">\(2 \pi \times \frac{ 2}{1024}\)</span> 和
<span class="math notranslate nohighlight">\(2 \pi \times \frac{3}{1024}\)</span> 也类似。
显然，一些重要的东西会因为精度太低而被丢弃。
因此，尽管输入数据的精度较低，但信号处理长期以来一直使用单精度浮点数来完成。
一种方式理解傅里叶变换的方式是把它视为波的和。
下面是一个具有两个波频率的 1024 点 FFT，实部显示为蓝色，虚部显示为紫色。</p>
<p>上述 FFT 代码的 unum 版本在附录 C.17 中。
与本章前面的示例不同，此计算涉及数万次计算，并为比较 unum 和 float
的存储效率提供了最佳基础。 首先，对答案进行健全性检查；
这是相同变换的图，但使用 unum 算术并跟踪位和个数的搬运量。</p>
<div class="figure align-default" id="id32">
<img alt="_images/image-20230628102943094.png" src="_images/image-20230628102943094.png" />
<p class="caption"><span class="caption-number">Fig. 236 </span><span class="caption-text">image-20230628102943094</span><a class="headerlink" href="#id32" title="Permalink to this image">¶</a></p>
</div>
<p>两个结果之间的差异远低于图表的分辨率； unum
版本中的边界宽度太小，在图中也看不到。 检查数据移动成本：</p>
<div class="figure align-default" id="id33">
<img alt="_images/image-20230628103711078.png" src="_images/image-20230628103711078.png" />
<p class="caption"><span class="caption-number">Fig. 237 </span><span class="caption-text">image-20230628103711078</span><a class="headerlink" href="#id33" title="Permalink to this image">¶</a></p>
</div>
<p>utag 只有 6 位长，最大 unum 大小为 25 位。 由于 FFT
中算术的许多或大部分结果都是
ubound，并且旋转因子引入了无理数，通常需要使用最大小数位数表示，因此您可能会以为每个数字的平均位数非常接近
ubound 有两个元素，即 51 位。 为什么每 ubound 使用的比特数少于 23
位？效果如此之好是怎么做到的。</p>
<p>大部分答案可以在 FFT
内核运算中找到，由于其数据流图的形状，有时被称为“蝶形运算”：</p>
<div class="figure align-default" id="id34">
<img alt="_images/image-20230628105422834.png" src="_images/image-20230628105422834.png" />
<p class="caption"><span class="caption-number">Fig. 238 </span><span class="caption-text">image-20230628105422834</span><a class="headerlink" href="#id34" title="Permalink to this image">¶</a></p>
</div>
<blockquote>
<div><center><p>“w”是“旋转因子”，对应于上述算法中的旋转</p>
</center></div></blockquote>
<p>假设 u 和 v 是 12 位整数，即 0 到 4095 之间的数字。在最上面的情况下，将
u 和 v 相加可能会导致整数溢出。 我们可能需要 13 位来准确表示总和。
在下面的情况下，我们从 u 中减去 v，这通常需要少于 12 位。
该算法实际上要求位级的灵活精度，因为固定大小的分数在一种情况下可能会出现舍入错误，而在另一种情况下可能会浪费存储空间。
unum 环境根据需要提升和降低精度，但仅此而已。</p>
<p>FFT 提供了尝试 <strong>smartunify</strong> 的机会。 原型 unum 函数使用目标比率为 1 的
<strong>smartunify</strong>，因此只有在增加每比特信息时才unify 含2部分的ubound。
如果不用它，前面 1K FFT 的成本将上升到每个值 28.5 位。
如图所示，边界的轻微放松对最终结果的影响可以忽略不计；
此外，unify后的答案仍然是有正确边界的，而不是像浮点计算那样的属于猜测。</p>
<p>石油和天然气行业长期以来一直依赖 32 位浮点来保护 12 位数据的中间计算，但
unum 方法提供了可证明正确的有界结果，每个数字的位数更少。
一旦用户对最终结果需要多少精度做出让步，根据应用的不同，还可以实现更大的节约。</p>
</div>
</div>
</div>
<div class="side-doc-outline">
<div class="side-doc-outline--content">
<div class="localtoc">
<p class="caption">
<span class="caption-text">Table Of Contents</span>
</p>
<ul>
<li><a class="reference internal" href="#">14 试运行：Unums 面临计算挑战</a><ul>
<li><a class="reference internal" href="#ii">14.1 浮点数II：卡汉之怒</a></li>
<li><a class="reference internal" href="#rump">14.2 Rump的烦恼事</a></li>
<li><a class="reference internal" href="#id1">14.3 二次公式</a></li>
<li><a class="reference internal" href="#id2">14.4 贝利的数字噩梦</a></li>
<li><a class="reference internal" href="#unum">14.5 使用unum做快速傅里叶变换</a></li>
</ul>
</li>
</ul>
</div>
</div>
</div>
<div class="clearer"></div>
</div><div class="pagenation">
<a id="button-prev" href="13_fused_operations.html" class="mdl-button mdl-js-button mdl-js-ripple-effect mdl-button--colored" role="botton" accesskey="P">
<i class="pagenation-arrow-L fas fa-arrow-left fa-lg"></i>
<div class="pagenation-text">
<span class="pagenation-direction">Previous</span>
<div>13 融合操作（一次性表达式）</div>
</div>
</a>
<a id="button-next" href="part1_summary.html" class="mdl-button mdl-js-button mdl-js-ripple-effect mdl-button--colored" role="botton" accesskey="N">
<i class="pagenation-arrow-R fas fa-arrow-right fa-lg"></i>
<div class="pagenation-text">
<span class="pagenation-direction">Next</span>
<div>小结</div>
</div>
</a>
</div>
</main>
</div>
</body>
</html>