\documentclass[12pt]{article} \usepackage{verbatim,moreverb,graphicx,amsmath} \def\listingoffset{1em} \def\listinglabel#1{\llap{\tiny\it\the#1}\hskip\listingoffset\relax} \def\mylisting#1{{\fontsize{10}{11}\selectfont\listinginput[1]{1}{#1}}} \def\myoutput#1{{\fontsize{10}{11}\selectfont\verbatimtabinput{#1}}} \newtheorem{problem}{Problem} \newenvironment{solution}[1][Solution]{\noindent{\it #1.} }{\ \rule {0.5em}{0.5em}\medskip\par} \newcommand{\artanh}{\mathop{\mathrm{artanh}}} \newcommand{\sinc}{\mathop{\mathrm{sinc}}} \def\words#1{\quad\hbox{#1}\quad} \def\wwords#1{\qquad\hbox{#1}\qquad} \def\Z{{\bf Z}} \def\R{{\bf R}} \def\P{{\bf P}} \def\x{{\tt x}} \def\ve{\varepsilon} \begin{document} \section{Homework 1} %\subsection{Chapter 2 Section 8} \begin{problem}[2.7.6] Convert $1/12$ to binary. \end{problem} \begin{solution} The doubling approach results in $$ \begin{matrix} 1/12&1/6&1/3&2/3&4/3&2/3&4/3&2/3& \ldots\\ . &0 &0 &0 &1 &0 &1 &0 & \ldots\\ \end{matrix} $$ from which we may infer that $$ {1\over 12} = 0.0001010\ldots = 0.000\overline{10}. $$ To verify the above formula note that \begin{align*} 0.000\overline{10}&=2^{-3}\cdot 0.\overline{10} =2^{-3}\sum_{k=1}^\infty 2\cdot2^{-2k}\\ &=2^{-2}\sum_{k=1}^\infty {1\over 4^k} =2^{-2}{1/4\over 1-1/4} =2^{-2}{1\over 3}={1\over 12}. \end{align*} \end{solution} \begin{problem}[2.7.7] Convert $1/16$ to binary. \end{problem} \begin{solution} Again use the doubling approach to obtain $$ \begin{matrix} 1/16&1/8&1/4&1/2&1&0&0&\ldots\\ . &0 &0 &0 &1&0&0&\ldots\\ \end{matrix} $$ which implies $$ {1\over 16}=0.0001. $$ \end{solution} \begin{problem}[2.7.11] Convert $0.111\,111$ to decimal. \end{problem} \begin{solution} Converting each group of three digits to base eight first yields $$ 0.111\,111={7\over 8}+{7\over 8^2}={56+7\over 64}=0.984375. $$ \end{solution} \begin{problem}[2.7.12] Convert $0.100\,001$ to decimal. \end{problem} Converting each group of three digits to base eight first yields $$ 0.100\,001={4\over 8}+{1\over 8^2}={32+1\over 64}=0.515625. $$ \begin{problem}[2.6.2] Examine the propagation of roundoff through division. \end{problem} \begin{solution} Suppose computations are performed using floating-point arithmetic accurate to $n$ significant digits. Let $x_1,x_2\in\R$ with $x_2\ne0$ and denote by $\x_1$ and $\x_2$ respectively the nearest floating-point approximations. Since by definition the relative errors in the approximations are bounded by $5\times 10^{-n}$, setting $$ \ve_1= {\x_1-x_1\over x_1} \wwords{and} \ve_2={\x_2-x_2\over x_2} $$ yields $\ve_1,\ve_2\in E$ where $E=[-5\times 10^{-n},5\times 10^{-n}]$. A simple rearrangement of these expressions yields that $$ \x_1=(1+\ve_1)x_1 \wwords{and} \x_2=(1+\ve_2)x_2. $$ Since $x_2\ne0$ it follows that $\x_2\ne 0$. Now, let $\x$ be the nearest floating-point approximation to the exact quotient $x=\x_1/\x_2$. Thus, $$ \x=(1+\ve)x\wwords{for some}\ve\in E. $$ We now combine the above to show how all the rounding errors have propagated into our estimate $\x$ of the exact quotient $x_1/x_2$. To do this compute the relative error as \begin{align*} {\x-x_1/x_2\over x_1/x_2} &={(1+\ve)x-x_1/x_2\over x_1/x_2} ={(1+\ve)\x_1/\x_2-x_1/x_2\over x_1/x_2}\\ &={(1+\ve)[(1+\ve_1)x_1]/[(1+\ve_2)x_2]-x_1/x_2\over x_1/x_2}\\ &={(1+\ve)(1+\ve_1)\over (1+\ve_2)}-1. \end{align*} In other words, $$ \x= {(1+\ve)(1+\ve_1)\over (1+\ve_2)}\cdot {x_1\over x_2}. $$ If we assume that all terms involving products of $\ve$, $\ve_1$ and $\ve_2$ can be neglected, the above expression can be simplified as \begin{align*} {(1+\ve)(1+\ve_1)\over (1+\ve_2)} &=(1+\ve+\ve_1+\ve\ve_1)\sum_{k=0}^\infty {-1\choose k} (\ve_2)^k\\ &\approx (1+\ve+\ve_1)(1-\ve_2) \approx (1+\ve+\ve_1-\ve_2). \end{align*} Therefore, the propagated relative error is approximately $\ve+\ve_1-\ve_2$ and is guaranteed to satisfy the bound $$ |\ve+\ve_1-\ve_2|\le 3\cdot 5\times 10^{-n}= 1.5\times 10^{-n+1}. $$ \end{solution} \begin{problem}[2.8.1] Derive the corresponding results for division relating to the frequency distribution of the mantissa. \end{problem} \begin{solution} We example the way division transforms various distributions. Let $x$ come from the probability density $f(x)$, let $y$ come from $g(y)$ and let the quotient $z$ have the density $h(z)$. Further denote the cumulative distributions by $$ F(x)=\int_{1/b}^x f(x)dx,\quad G(y)=\int_{1/b}^y g(y)dy\words{and} H(z)=\int_{1/b}^z h(z)dz. $$ In order to understand how $H(z)$ depends on the distributions $f(x)$ and $g(y)$ we first note the mantissa of the quotient corresponding to $x/y$ is given by $$ m(x,y)=\begin{cases} x/y&x/y<1\cr x/(yb)&x/y \ge 1.\cr \end{cases} $$ where second case occurs when the exponent has to be shifted by one so $z$ is still of the form{\tt\ 0.DDDDDDDDD }where the{\tt\ D}'s represent the base-$b$ digits the first of which is non-zero. We interpret $H(z)=\P\{ m(x,y)\le z \}$ as the probability measure of the set of points for which the mantissa of the $x/y$ is less than or equal $z$. Graphically the condition $m(x,y)\le z$ may be depicted by $$ \includegraphics{gp281.pdf} $$ We now integrate over this set to obtain the cumulative distribution function \begin{align*} H(z)&=\int_{1/b}^z\int_{x/z}^1 f(x)g(y)dy\,dx +\int_{1/b}^z\int_{1/b}^x f(x)g(y)dy\,dx\\ &\qquad\qquad+\int_z^1\int_{x/(zb)}^x f(x)g(y)dy\,dx\\ &=\int_{1/b}^z f(x)\Big[G(1)-G\Big({x\over z}\Big)\Big]dx +\int_{1/b}^z f(x)\Big[G(x)-G\Big({1\over b}\Big)\Big]dx\\ &\qquad\qquad+\int_z^1 f(x)\Big[G(x)-G\Big({x\over zb}\Big)\Big]dx. \end{align*} Then differentiating gives the density \begin{align*} h(z)&=f(z)\Big[G(1)-G\Big({z\over z}\Big) +G(z)-G\Big({1\over b}\Big)-G(z)+G\Big({z\over zb}\Big) \Big]\\ &+\int_{1/b}^z f(x)g\Big({x\over z}\Big){x\over z^2}dx +{1\over b}\int_z^1 f(x)g\Big({x\over zb}\Big){x\over z^2}dx. \end{align*} Suppose the divisor has the reciprocal distribution; that is, $$ g(y)={1\over y\log b}. $$ Substituting this into the formula for $h(z)$ yields \begin{align*} h(z)&= \int_{1/b}^z f(x)\Big({z\over x\log b}\Big){x\over z^2}dx +{1\over b}\int_z^1 f(x)\Big({zb\over x\log b}\Big){x\over z^2}dx\\ &=\int_{1/b}^z f(x)\Big({1\over z\log b}\Big)dx +\int_z^1 f(x)\Big({1\over z\log b}\Big)dx ={1\over z\log b}. \end{align*} Thus, if one of the divisor of a quotient comes from the reciprocal distribution, then regardless of the distribution of the dividend, the quotient has the reciprocal distribution. Note for products treated in the text there is a symmetry which implies one does not need consider separately the case when $f(x)$ is given by the reciprocal distribution. However, such is not the case with quotients, so we alternatively suppose $$ f(x)={1\over x\log b}. $$ Now, substituting this into the formula for $h(z)$ yields \begin{align*} h(z) &=\int_{1/b}^z {1\over x\log b}g\Big({x\over z}\Big){x\over z^2}dx +{1\over b}\int_z^1 {1\over x\log b}g\Big({x\over zb}\Big){x\over z^2}dx\\ &={1\over z\log b}\bigg\{ \int_{1/b}^z g\Big({x\over z}\Big){dx\over z} +\int_z^1 g\Big({x\over zb}\Big){dx\over zb}\bigg\}\\ &={1\over z\log b}\bigg\{ \int_{1/(zb)}^1 g(u)du +\int_{1/b}^{1/(zb)} g(v)dv\bigg\} ={1\over z\log b}. \end{align*} Here we have made the substitutions $u=x/z$ and $v=x/(zb)$. Having treated both cases, we see that provided either $f(x)$ or $g(y)$ have the reciprocal distribution then $h(z)$ will also have the reciprocal distribution. Thus, the persistence of the reciprocal distribution is also established in the case of division. Next we show how the reciprocal distribution can arise. Let $r(z)$ be the reciprocal distribution and define $$ D\{h\}= \max\Big\{ \Big|{h(z)-r(z)\over r(z)}\Big| : {1\over b}\le z\le 1 \Big\}. $$ Since $$ r(z)= \int_{1/b}^z f(x)r\Big({x\over z}\Big){x\over z^2}dx +{1\over b}\int_z^1 f(x)r\Big({x\over zb}\Big){x\over z^2}dx, $$ it follows that \begin{align*} {h(z)-r(z)\over r(z)} &= \int_{1/b}^z f(x)\Big\{{g(x/z)-r(x/z)\over r(z)}\Big\} {x\over z^2}dx\\ &+ {1\over b}\int_z^1 f(x) \Big\{ {g[x/(zb)]-r[x/(zb)]\over r(z)} \Big\} {x\over z^2}dx. \end{align*} But $$ {z^2\over x} r(z)= {z\over x\log b} = r\Big({x\over z}\Big) \wwords{and} {bz^2\over x} r(z)= {bz\over x\log b} = r\Big({x\over bz}\Big) $$ imply \begin{align*} {h(z)-r(z)\over r(z)} &= \int_{1/b}^z f(x)\Big\{{g(x/z)-r(x/z)\over r(x/z)}\Big\} dx\\ &+ \int_z^1 f(x) \Big\{ {g[x/(zb)]-r[x/(zb)]\over r[x/(zb)]} \Big\} dx. \end{align*} Since $f(x)\ge 0$ in the two intervals $$ \Big|{h(z)-r(z)\over r(z)}\Big| \le \int_{1/b}^z f(x) D\{g\} dx + \int_z^1 f(x) D\{g\} dx \le D\{g\}. $$ for all z. Consequently $D\{h\}\le D\{g\}$. By a similar argument $D\{h\}\le D\{f\}$. In particular, it is impossible for the distribution of mantissas of the quotient to get farther from the reciprocal distribution than either the distribution of the dividend or the divisor. \end{solution} \begin{problem}[3.3.1] For large $x$ rearrange $$ {1\over 1+x}-{1\over x} $$ for evaluation with respect to $x$. \end{problem} \begin{solution} Simply by finding a common denominator as $$ {1\over 1+x}-{1\over x}={x-(1+x)\over x(1+x)} ={-1\over x(1+x)}. $$ Note after simplifying that the difference of two nearly equal numbers has been eliminated from the expression so there is no longer any loss of precision due to cancellation. \end{solution} \begin{problem}[3.3.4] For large $x$ rearrange $$ \root 3\of{x+1} - \root 3\of{x} $$ for evaluation with respect to $x$. \end{problem} \begin{solution} Since $a^3-b^3=(a-b)(a^2+ab+b^2)$ then simplify by completing the difference of cubes as \begin{align*} \root 3\of{x+1} - \root 3\of{x} &= \big((x+1)^{1/3} - x^{1/3}\big) \cdot{ (x+1)^{2/3} + (x+1)^{1/3}x^{1/3}+x^{2/3}\over (x+1)^{2/3} + (x+1)^{1/3}x^{1/3}+x^{2/3}}\\ &= {1\over (x+1)^{2/3} + (x+1)^{1/3}x^{1/3}+x^{2/3}}. \end{align*} After simplifying the difference of two nearly equal numbers has been eliminated so there is no longer any loss of precision due to cancellation. \end{solution} \begin{problem}[3.3.5] For small $\ve$ rearrange $$ \cos(x+\ve)-\cos x $$ for evaluation with respect to $x$. \end{problem} \begin{solution} Using the angle addition and subtraction formulas \begin{align*} % \sin(a+b)&=\sin a\cos b+\cos a\sin b\\ \cos(a+b)&=\cos a\cos b-\sin a\sin b\\ \cos(a-b)&=\cos a\cos b+\sin a\sin b \end{align*} write \begin{align*} \cos(x+\ve)&= \cos\big((x+\ve/2)+\ve/2\big)\\ &=\cos(x+\ve/2)\cos(\ve/2) -\sin(x+\ve/2)\sin(\ve/2)\\ \cos(x)&= \cos\big((x+\ve/2)-\ve/2\big)\\ &=\cos(x+\ve/2)\cos(\ve/2) +\sin(x+\ve/2)\sin(\ve/2). \end{align*} Consequently \begin{align*} \cos(x+\ve)-\cos x &=-2 \sin(x+\ve/2)\sin(\ve/2). \end{align*} Note the subtraction of two nearly equal numbers when $\ve$ is small has been eliminated; however, if $\sin x\approx 0$ then loss of precision may still occur in the evaluation of $\sin(x+\ve/2)$. In this case $x\approx k\pi$ for some $k\in\Z$ and there is no way that I know to avoid this other way of losing precision. \end{solution} \begin{problem}[3.3.6] For large $N$ rearrange $$ \int_N^{N+1} {dx\over x} = \ln(N+1)-\ln N $$ for evaluation with respect to $N$. \end{problem} \begin{solution} Recall the definition of the inverse hyperbolic tangent $$ \artanh x = {1\over 2}\ln\Big({1+x\over 1-x}\Big) $$ which is available as the standard subroutine{\tt\ atanh }for evaluation on the computer. Setting $x=1/(2N+1)$ yields $$ {1+x\over 1-x} ={1+1/(2N+1)\over 1-1/(2N+1)} ={2N+2\over 2N}={N+1\over N}. $$ Therefore $$ \int_N^{N+1}{dx\over x}=\ln(N+1)-\ln N =\ln\Big({N+1\over N}\Big) =2\artanh\Big({1\over 2N+1}\Big). $$ Note that a direct evaluation of the expression $$ \ln\Big(1+{1\over N}\Big) %\approx {\tt log}\Big(1+{1\over N}\Big) $$ on the computer using the built-in{\tt\ log }function for $N$ large also leads to loss of precision. This is why we used the inverse hyperbolic tangent. The computer may also have a function called{\tt\ log1p }defined as ${\tt log1p}(x)=\log(1+x)$ which could be used instead. In this case one could simply evaluate $$ \int_N^{N+1}{dx\over x}=\ln(N+1)-\ln N \approx {\tt log1p}\Big({1\over N}\Big) $$ without any loss of precision. \end{solution} \begin{problem}[3.4.1] For $\ve$ small compute $$ {1-\cos\ve\over \ve^2}. $$ \end{problem} \begin{solution} It's tempting to write $$ {1-\cos\ve\over \ve^2} ={1-\cos\ve\over \ve^2}\cdot {1+\cos\ve\over 1+\cos\ve} =\Big({\sin\ve\over \ve}\Big)^2 \Big({1\over1+\cos\ve}\Big) = {\sinc^2\ve\over1+\cos\ve}. $$ This removes the loss of precision, but unfortunately{\tt sinc }is not part of the standard mathematical library available for evaluation on the computer. One could use an approximation such as $$ \sinc\ve=\begin{cases} (\sin \ve)/\ve & \mbox{for $|\ve|>10^{-15}$}\cr 0 & \mbox{otherwise}, \end{cases} $$ but since this section is about using series we choose a different approach. Since $$ \cos \ve = \sum_{k=0}^\infty {(-1)^k\over (2k)!}\ve^{2k} $$ then $$ 1-\cos \ve = -\sum_{k=1}^\infty {(-1)^k\over (2k)!}\ve^{2k} = \sum_{k=0}^\infty {(-1)^k\over (2k+2)!}\ve^{2k+2} $$ and consequently $$ {1-\cos\ve\over\ve^2} = \sum_{k=0}^\infty {(-1)^k\over (2k+2)!}\ve^{2k} \approx {\frac{1}{2}}-{\frac{1}{24}}{x}^{2}+{\frac{1}{720}}{x}^{4}-{\frac{1}{ 40320}}{x}^{6}+\cdots. $$ Supposing $\ve$ is small one may now truncate at only a few terms $N$ such that $$ {2\over (2N+2)!}\ve^{2N} < 10^{-15} $$ to obtain an approximation good to within the accuracy of double-precision floating-point arithmetic. \end{solution} \begin{problem}[3.4.4] For $\ve$ small compute $$ \sqrt{{e^{2\ve}-1\over e^\ve-1}}. $$ \end{problem} \begin{solution} Again it is tempting to make the exact algebraic simplifications $$ {e^{2\ve}-1\over e^\ve-1}= {(e^\ve-1)(e^\ve+1)\over e^\ve-1}= e^\ve+1 $$ and the compute $$ \sqrt{{e^{2\ve}-1\over e^\ve-1}} =\sqrt{e^\ve+1}. $$ This results in a exact expression which avoids any loss of precision due to subtraction of nearly equal numbers. However, this expression still involves the numerical computation of{\tt\ exp }and{\tt\ sqrt }which are unnecessarily time consuming when $\ve$ is known to be small. Rather than expanding the original expression in terms of series we proceed as follows. $$ e^\ve+1 =\Big(\sum_{k=0}^\infty {1\over k!}\ve^k\Big) + 1 =2+\sum_{k=1}^\infty {1\over k!}\ve^k =2(1+\delta) $$ where $$ \delta={1\over 2}\sum_{k=1}^\infty {1\over k!}\ve^k. $$ By the binomial theorem \begin{align*} \sqrt{e^\ve+1}&=\sqrt{2}\cdot\sqrt{1+\delta} =\sqrt{2}\sum_{l=0}^\infty {1/2\choose l}\delta^l\\ &=\sqrt{2}\sum_{l=0}^\infty {1/2\choose l} \Big( {1\over 2}\sum_{k=1}^\infty {1\over k!}\ve^k \Big)^l\\ &=\sqrt 2 \Big(1+{1\over 4}x+{\frac {3}{32}} {x}^{2} +{\frac {7}{384}}{x}^{3} +{\frac {3}{2048}}{x}^{4}+\cdots\Big). \end{align*} \end{solution} \begin{problem}[3.4.5] For $\ve$ small find an expression to compute $$ {\ve-\sin\ve\over \ve-\tan\ve}. $$ \end{problem} \begin{solution} First note that \begin{align*} {\ve-\sin\ve\over \ve-\tan\ve} ={\ve \cos\ve -\sin\ve\cos\ve \over \ve \cos\ve-\sin\ve} ={\ve \cos\ve -2^{-1}\sin(2\ve) \over \ve \cos\ve-\sin\ve}. \end{align*} Now \begin{align*} \ve\cos\ve-\sin\ve &= \ve\sum_{k=0}^\infty {(-1)^k\over (2k)!}\ve^{2k} -\sum_{k=0}^\infty {(-1)^k\over (2k+1)!}\ve^{2k+1}\\ &=\ve\sum_{k=0}^\infty (-1)^k \Big({1\over (2k)!}-{1\over (2k+1)!}\Big)\ve^{2k}\\ &=\ve\sum_{k=1}^\infty { (-1)^k 2k\over (2k+1)!}\ve^{2k} \end{align*} and similarly \begin{align*} \ve\cos\ve-2^{-1}\sin(2\ve) &= \ve\sum_{j=0}^\infty {(-1)^j\over (2j)!}\ve^{2j} -{1\over 2}\sum_{j=0}^\infty {(-1)^j\over (2j+1)!}(2\ve)^{2j+1}\\ &=\ve\sum_{j=0}^\infty (-1)^j \Big({1\over (2j)!}-{2^{2j}\over (2j+1)!}\Big)\ve^{2j}\\ &= \ve\sum_{j=1}^\infty { (-1)^j (2j+1-2^{2j})\over (2j+1)!}\ve^{2j}. \end{align*} By the binomial theorem \begin{align*} {\ve^3\over \ve\cos\ve-\sin\ve} &= {1 \over \sum_{k=1}^\infty { (-1)^k 2k\over (2k+1)!}\ve^{2k-2}}\\ &={1 \over -{1\over 3}+\sum_{k=2}^\infty { (-1)^k 2k\over (2k+1)!}\ve^{2k-2}} = {-3\over 1+\delta} =-3 \sum_{l=0}^\infty {-1\choose l} \delta^l \end{align*} where $$ \delta=-3 \sum_{k=2}^\infty { (-1)^k 2k\over (2k+1)!}\ve^{2k-2}. $$ It follows that \begin{align*} &{\ve\cos\ve-2^{-1}\sin(2\ve) \over \ve\cos\ve-\sin\ve} = -3 \sum_{j=1}^\infty { (-1)^j (2j+1-2^{2j})\over (2j+1)!}\ve^{2j-2} \sum_{l=0}^\infty {-1\choose l} \delta^l\\ &= -3 \sum_{j=1}^\infty { (-1)^j (2j+1-2^{2j})\over (2j+1)!}\ve^{2j-2} \sum_{l=0}^\infty {-1\choose l} \Big(-3 \sum_{k=2}^\infty { (-1)^k 2k\over (2k+1)!}\ve^{2k-2}\Big)^l\\ &=-{\frac {1}{2}}+{\frac {9}{40}}{{\it \ve}}^{2}-{\frac {27}{2800}}{{ \it \ve}}^{4}+{\frac {27}{112000}}{{\it \ve}}^{6}-\cdots. \end{align*} \end{solution} \begin{problem}[3.4.6] For $x$ large find an expression to compute $$ (x+1)^{1/n}-x^{1/n}. $$ \end{problem} \begin{solution} Factoring and applying the binomial theorem yields \begin{align*} (x+1)^{1/n}-x^{1/n} &= x^{1/n}\big( (1+1/x)^{1/n}-1\big)\\ &= x^{1/n}\bigg\{\Big( \sum_{l=0}^\infty {1/n\choose l} (1/x)^l\Big) -1\bigg\}\\ &= x^{1/n}\sum_{l=1}^\infty {1/n\choose l} (1/x)^l\\ &= x^{1/n}\Big( {1\over nx}+{n-1\over 2(nx)^2} +{(n-1)(2n-1)\over 6(nx)^3}\\ &\qquad\qquad+{(n-1)(2n-1)(3n-1)\over 24(nx)^4}+\cdots \Big). \end{align*} \end{solution} \begin{problem}[4.7.1] Apply Newton's method using $x_0=1/2$ to solve $$ xe^x-1=0. $$ \end{problem} \begin{solution} After computing $$ f(x)=xe^x-1\wwords{and} f'(x)=(1+x)e^x, $$ I wrote a program in $C$ that looked like \mylisting{p471.c} The output was \myoutput{p471.out} This shows $x\approx 0.56714329040978$. \end{solution} \begin{problem}[4.7.2] Apply Newton's method using $x_0=1$ to solve $$ \arctan(x)-1=0. $$ \end{problem} \begin{solution} After computing $$ f(x)=\arctan(x)-1\wwords{and} f'(x)={1\over 1+x^2}, $$ I wrote a program in $C$ that looked like \mylisting{p472.c} The output was \myoutput{p472.out} This shows $x\approx 1.5574077246549$. \end{solution} \begin{problem}[4.7.3] Apply Newton's method using $x_0=10$ to solve $$ \ln(x) -3 = 0. $$ \end{problem} \begin{solution} After computing $$ f(x)=\log(x)-3\wwords{and} f'(x)={1\over x}, $$ I wrote a program in $C$ that looked like \mylisting{p473.c} The output was \myoutput{p473.out} This shows $x\approx 20.085536923188$. \end{solution} \end{document} The graph was this: \includegraphics[width=5.5in]{poly.eps} \end{document}