| 1 | <!DOCTYPE html PUBLIC "-//W3C//DTD HTML 4.01 Transitional//EN" "http://www.w3.org/TR/html4/loose.dtd">
|
|---|
| 2 | <html>
|
|---|
| 3 | <!-- Copyright (C) 2022 Richard Stallman and Free Software Foundation, Inc.
|
|---|
| 4 |
|
|---|
| 5 | (The work of Trevis Rothwell and Nelson Beebe has been assigned or
|
|---|
| 6 | licensed to the FSF.)
|
|---|
| 7 |
|
|---|
| 8 | Permission is granted to copy, distribute and/or modify this document
|
|---|
| 9 | under the terms of the GNU Free Documentation License, Version 1.3 or
|
|---|
| 10 | any later version published by the Free Software Foundation; with the
|
|---|
| 11 | Invariant Sections being "GNU General Public License," with the
|
|---|
| 12 | Front-Cover Texts being "A GNU Manual," and with the Back-Cover
|
|---|
| 13 | Texts as in (a) below. A copy of the license is included in the
|
|---|
| 14 | section entitled "GNU Free Documentation License."
|
|---|
| 15 |
|
|---|
| 16 | (a) The FSF's Back-Cover Text is: "You have the freedom to copy and
|
|---|
| 17 | modify this GNU manual. Buying copies from the FSF supports it in
|
|---|
| 18 | developing GNU and promoting software freedom." -->
|
|---|
| 19 | <!-- Created by GNU Texinfo 6.7, http://www.gnu.org/software/texinfo/ -->
|
|---|
| 20 | <head>
|
|---|
| 21 | <meta http-equiv="Content-Type" content="text/html; charset=utf-8">
|
|---|
| 22 | <title>Error Recovery (GNU C Language Manual)</title>
|
|---|
| 23 |
|
|---|
| 24 | <meta name="description" content="Error Recovery (GNU C Language Manual)">
|
|---|
| 25 | <meta name="keywords" content="Error Recovery (GNU C Language Manual)">
|
|---|
| 26 | <meta name="resource-type" content="document">
|
|---|
| 27 | <meta name="distribution" content="global">
|
|---|
| 28 | <meta name="Generator" content="makeinfo">
|
|---|
| 29 | <link href="index.html" rel="start" title="Top">
|
|---|
| 30 | <link href="Symbol-Index.html" rel="index" title="Symbol Index">
|
|---|
| 31 | <link href="index.html#SEC_Contents" rel="contents" title="Table of Contents">
|
|---|
| 32 | <link href="Floating-Point-in-Depth.html" rel="up" title="Floating Point in Depth">
|
|---|
| 33 | <link href="Exact-Floating-Constants.html" rel="next" title="Exact Floating Constants">
|
|---|
| 34 | <link href="Fused-Multiply_002dAdd.html" rel="prev" title="Fused Multiply-Add">
|
|---|
| 35 | <style type="text/css">
|
|---|
| 36 | <!--
|
|---|
| 37 | a.summary-letter {text-decoration: none}
|
|---|
| 38 | blockquote.indentedblock {margin-right: 0em}
|
|---|
| 39 | div.display {margin-left: 3.2em}
|
|---|
| 40 | div.example {margin-left: 3.2em}
|
|---|
| 41 | div.lisp {margin-left: 3.2em}
|
|---|
| 42 | kbd {font-style: oblique}
|
|---|
| 43 | pre.display {font-family: inherit}
|
|---|
| 44 | pre.format {font-family: inherit}
|
|---|
| 45 | pre.menu-comment {font-family: serif}
|
|---|
| 46 | pre.menu-preformatted {font-family: serif}
|
|---|
| 47 | span.nolinebreak {white-space: nowrap}
|
|---|
| 48 | span.roman {font-family: initial; font-weight: normal}
|
|---|
| 49 | span.sansserif {font-family: sans-serif; font-weight: normal}
|
|---|
| 50 | ul.no-bullet {list-style: none}
|
|---|
| 51 | -->
|
|---|
| 52 | </style>
|
|---|
| 53 |
|
|---|
| 54 |
|
|---|
| 55 | </head>
|
|---|
| 56 |
|
|---|
| 57 | <body lang="en">
|
|---|
| 58 | <span id="Error-Recovery"></span><div class="header">
|
|---|
| 59 | <p>
|
|---|
| 60 | Next: <a href="Exact-Floating-Constants.html" accesskey="n" rel="next">Exact Floating Constants</a>, Previous: <a href="Fused-Multiply_002dAdd.html" accesskey="p" rel="prev">Fused Multiply-Add</a>, Up: <a href="Floating-Point-in-Depth.html" accesskey="u" rel="up">Floating Point in Depth</a> [<a href="index.html#SEC_Contents" title="Table of contents" rel="contents">Contents</a>][<a href="Symbol-Index.html" title="Index" rel="index">Index</a>]</p>
|
|---|
| 61 | </div>
|
|---|
| 62 | <hr>
|
|---|
| 63 | <span id="Error-Recovery-1"></span><h3 class="section">28.11 Error Recovery</h3>
|
|---|
| 64 | <span id="index-error-recovery-_0028floating-point_0029"></span>
|
|---|
| 65 | <span id="index-floating_002dpoint-error-recovery"></span>
|
|---|
| 66 |
|
|---|
| 67 | <p>When two numbers are combined by one of the four basic
|
|---|
| 68 | operations, the result often requires rounding to storage
|
|---|
| 69 | precision. For accurate computation, one would like to be able
|
|---|
| 70 | to recover that rounding error. With historical floating-point
|
|---|
| 71 | designs, it was difficult to do so portably, but now that IEEE
|
|---|
| 72 | 754 arithmetic is almost universal, the job is much easier.
|
|---|
| 73 | </p>
|
|---|
| 74 | <p>For addition with the default <em>round-to-nearest</em> rounding
|
|---|
| 75 | mode, we can determine the error in a sum like this:
|
|---|
| 76 | </p>
|
|---|
| 77 | <div class="example">
|
|---|
| 78 | <pre class="example">volatile double err, sum, tmp, x, y;
|
|---|
| 79 |
|
|---|
| 80 | if (fabs (x) >= fabs (y))
|
|---|
| 81 | {
|
|---|
| 82 | sum = x + y;
|
|---|
| 83 | tmp = sum - x;
|
|---|
| 84 | err = y - tmp;
|
|---|
| 85 | }
|
|---|
| 86 | else /* fabs (x) < fabs (y) */
|
|---|
| 87 | {
|
|---|
| 88 | sum = x + y;
|
|---|
| 89 | tmp = sum - y;
|
|---|
| 90 | err = x - tmp;
|
|---|
| 91 | }
|
|---|
| 92 | </pre></div>
|
|---|
| 93 |
|
|---|
| 94 | <p><span id="index-twosum"></span>
|
|---|
| 95 | Now, <code>x + y</code> is <em>exactly</em> represented by <code>sum + err</code>.
|
|---|
| 96 | This basic operation, which has come to be called <em>twosum</em>
|
|---|
| 97 | in the numerical-analysis literature, is the first key to tracking,
|
|---|
| 98 | and accounting for, rounding error.
|
|---|
| 99 | </p>
|
|---|
| 100 | <p>To determine the error in subtraction, just swap the <code>+</code> and
|
|---|
| 101 | <code>-</code> operators.
|
|---|
| 102 | </p>
|
|---|
| 103 | <p>We used the <code>volatile</code> qualifier (see <a href="volatile.html">volatile</a>) in the
|
|---|
| 104 | declaration of the variables, which forces the compiler to store and
|
|---|
| 105 | retrieve them from memory, and prevents the compiler from optimizing
|
|---|
| 106 | <code>err = y - ((x + y) - x)</code> into <code>err = 0</code>.
|
|---|
| 107 | </p>
|
|---|
| 108 | <p>For multiplication, we can compute the rounding error without
|
|---|
| 109 | magnitude tests with the FMA operation (see <a href="Fused-Multiply_002dAdd.html">Fused Multiply-Add</a>),
|
|---|
| 110 | like this:
|
|---|
| 111 | </p>
|
|---|
| 112 | <div class="example">
|
|---|
| 113 | <pre class="example">volatile double err, prod, x, y;
|
|---|
| 114 | prod = x * y; /* <span class="roman">rounded product</span> */
|
|---|
| 115 | err = fma (x, y, -prod); /* <span class="roman">exact product <code>= <var>prod</var> + <var>err</var></code></span> */
|
|---|
| 116 | </pre></div>
|
|---|
| 117 |
|
|---|
| 118 | <p>For addition, subtraction, and multiplication, we can represent the
|
|---|
| 119 | exact result with the notional sum of two values. However, the exact
|
|---|
| 120 | result of division, remainder, or square root potentially requires an
|
|---|
| 121 | infinite number of digits, so we can at best approximate it.
|
|---|
| 122 | Nevertheless, we can compute an error term that is close to the true
|
|---|
| 123 | error: it is just that error value, rounded to machine precision.
|
|---|
| 124 | </p>
|
|---|
| 125 | <p>For division, you can approximate <code>x / y</code> with <code>quo + err</code>
|
|---|
| 126 | like this:
|
|---|
| 127 | </p>
|
|---|
| 128 | <div class="example">
|
|---|
| 129 | <pre class="example">volatile double err, quo, x, y;
|
|---|
| 130 | quo = x / y;
|
|---|
| 131 | err = fma (-quo, y, x) / y;
|
|---|
| 132 | </pre></div>
|
|---|
| 133 |
|
|---|
| 134 | <p>For square root, we can approximate <code>sqrt (x)</code> with <code>root +
|
|---|
| 135 | err</code> like this:
|
|---|
| 136 | </p>
|
|---|
| 137 | <div class="example">
|
|---|
| 138 | <pre class="example">volatile double err, root, x;
|
|---|
| 139 | root = sqrt (x);
|
|---|
| 140 | err = fma (-root, root, x) / (root + root);
|
|---|
| 141 | </pre></div>
|
|---|
| 142 |
|
|---|
| 143 | <p>With the reliable and predictable floating-point design provided
|
|---|
| 144 | by IEEE 754 arithmetic, we now have the tools we need to track
|
|---|
| 145 | errors in the five basic floating-point operations, and we can
|
|---|
| 146 | effectively simulate computing in twice working precision, which
|
|---|
| 147 | is sometimes sufficient to remove almost all traces of arithmetic
|
|---|
| 148 | errors.
|
|---|
| 149 | </p>
|
|---|
| 150 |
|
|---|
| 151 |
|
|---|
| 152 |
|
|---|
| 153 | <hr>
|
|---|
| 154 | <div class="header">
|
|---|
| 155 | <p>
|
|---|
| 156 | Next: <a href="Exact-Floating-Constants.html" accesskey="n" rel="next">Exact Floating Constants</a>, Previous: <a href="Fused-Multiply_002dAdd.html" accesskey="p" rel="prev">Fused Multiply-Add</a>, Up: <a href="Floating-Point-in-Depth.html" accesskey="u" rel="up">Floating Point in Depth</a> [<a href="index.html#SEC_Contents" title="Table of contents" rel="contents">Contents</a>][<a href="Symbol-Index.html" title="Index" rel="index">Index</a>]</p>
|
|---|
| 157 | </div>
|
|---|
| 158 |
|
|---|
| 159 |
|
|---|
| 160 |
|
|---|
| 161 | </body>
|
|---|
| 162 | </html>
|
|---|
