diff options
| -rw-r--r-- | .DS_Store | bin | 0 -> 6148 bytes | |||
| -rw-r--r-- | homeworks/hw-2.org | 28 | ||||
| -rw-r--r-- | homeworks/hw-2.pdf | bin | 68368 -> 68190 bytes | |||
| -rw-r--r-- | homeworks/hw-2.tex | 36 | ||||
| -rw-r--r-- | homeworks/hw-3.org | 244 | ||||
| -rw-r--r-- | homeworks/hw-3.pdf | bin | 0 -> 152997 bytes | |||
| -rw-r--r-- | homeworks/hw-3.tex | 250 | ||||
| -rw-r--r-- | notes/Oct-4.org | 22 | ||||
| -rw-r--r-- | notes/Oct-6.org | 13 | ||||
| -rw-r--r-- | notes/Sep-15.org | 8 | ||||
| -rw-r--r-- | notes/Sep-15.pdf | bin | 0 -> 123321 bytes | |||
| -rw-r--r-- | notes/Sep-15.tex | 88 | ||||
| -rw-r--r-- | src/approx/derivative.lisp | 10 | ||||
| -rw-r--r-- | src/approx/maceps.lisp | 6 | ||||
| -rw-r--r-- | src/approx/package.lisp | 3 | ||||
| -rw-r--r-- | tests/approx.lisp | 16 |
16 files changed, 670 insertions, 54 deletions
diff --git a/.DS_Store b/.DS_Store Binary files differnew file mode 100644 index 0000000..e458869 --- /dev/null +++ b/.DS_Store diff --git a/homeworks/hw-2.org b/homeworks/hw-2.org index 29502d7..f696ffb 100644 --- a/homeworks/hw-2.org +++ b/homeworks/hw-2.org @@ -20,35 +20,10 @@ Computing $\epsilon_{\text{mac}}$ for single precision numbers :domain-values domain-values))) #+END_SRC -#+RESULTS: -| a | h | err | -| 1.0 | 0.5 | 0.5 | -| 1.0 | 0.25 | 0.25 | -| 1.0 | 0.125 | 0.125 | -| 1.0 | 0.0625 | 0.0625 | -| 1.0 | 0.03125 | 0.03125 | -| 1.0 | 0.015625 | 0.015625 | -| 1.0 | 0.0078125 | 0.0078125 | -| 1.0 | 0.00390625 | 0.00390625 | -| 1.0 | 0.001953125 | 0.001953125 | -| 1.0 | 0.0009765625 | 0.0009765625 | -| 1.0 | 0.00048828125 | 0.00048828125 | -| 1.0 | 0.00024414063 | 0.00024414063 | -| 1.0 | 0.00012207031 | 0.00012207031 | -| 1.0 | 6.1035156e-05 | 6.1035156e-05 | -| 1.0 | 3.0517578e-05 | 3.0517578e-05 | -| 1.0 | 1.5258789e-05 | 1.5258789e-05 | -| 1.0 | 7.6293945e-06 | 7.6293945e-06 | -| 1.0 | 3.8146973e-06 | 3.8146973e-06 | -| 1.0 | 1.9073486e-06 | 1.9073486e-06 | -| 1.0 | 9.536743e-07 | 9.536743e-07 | -| 1.0 | 4.7683716e-07 | 4.7683716e-07 | -| 1.0 | 2.3841858e-07 | 2.3841858e-07 | -| 1.0 | 1.1920929e-07 | 1.1920929e-07 | - (with many rows truncated) | a | h | err | +| 1.0 | 1.0 | 1.0 | | 1.0 | 0.5 | 0.5 | | 1.0 | 0.25 | 0.25 | | 1.0 | 0.125 | 0.125 | @@ -76,6 +51,7 @@ Computing $\epsilon_{\text{mac}}$ for double precision numbers: (with many rows truncated) | a | h | err | +| 1.0d0 | 1.0d0 | 1.0d0 | | 1.0d0 | 0.5d0 | 0.5d0 | | 1.0d0 | 0.25d0 | 0.25d0 | | 1.0d0 | 0.125d0 | 0.125d0 | diff --git a/homeworks/hw-2.pdf b/homeworks/hw-2.pdf Binary files differindex b719796..2dc4d28 100644 --- a/homeworks/hw-2.pdf +++ b/homeworks/hw-2.pdf diff --git a/homeworks/hw-2.tex b/homeworks/hw-2.tex index e0aa172..da8d3f5 100644 --- a/homeworks/hw-2.tex +++ b/homeworks/hw-2.tex @@ -1,4 +1,4 @@ -% Created 2023-09-27 Wed 10:09 +% Created 2023-10-07 Sat 14:51 % Intended LaTeX compiler: pdflatex \documentclass[11pt]{article} \usepackage[utf8]{inputenc} @@ -29,7 +29,7 @@ \setlength\parindent{0pt} \section{Question One} -\label{sec:orga21c813} +\label{sec:org58b9af4} Computing \(\epsilon_{\text{mac}}\) for single precision numbers \begin{verbatim} @@ -49,6 +49,7 @@ Computing \(\epsilon_{\text{mac}}\) for single precision numbers \begin{center} \begin{tabular}{rrr} a & h & err\\[0pt] +1.0 & 1.0 & 1.0\\[0pt] 1.0 & 0.5 & 0.5\\[0pt] 1.0 & 0.25 & 0.25\\[0pt] 1.0 & 0.125 & 0.125\\[0pt] @@ -65,7 +66,7 @@ a & h & err\\[0pt] \(\epsilon_{\text{mac single precision}}\) \(\approx\) 1.192(10\textsuperscript{-7}) \section{Question Two} -\label{sec:org06c4a23} +\label{sec:org27557b4} Computing \(\epsilon_{\text{mac}}\) for double precision numbers: \begin{verbatim} @@ -81,6 +82,7 @@ Computing \(\epsilon_{\text{mac}}\) for double precision numbers: \begin{center} \begin{tabular}{rrr} a & h & err\\[0pt] +1.0d0 & 1.0d0 & 1.0d0\\[0pt] 1.0d0 & 0.5d0 & 0.5d0\\[0pt] 1.0d0 & 0.25d0 & 0.25d0\\[0pt] 1.0d0 & 0.125d0 & 0.125d0\\[0pt] @@ -102,7 +104,7 @@ a & h & err\\[0pt] Thus, \(\epsilon_{\text{mac double precision}}\) \(\approx\) 2.220 \(\cdot\) 10\textsuperscript{-16} \section{Question Three - |v|\textsubscript{2}} -\label{sec:orgf181ba7} +\label{sec:org59c6c10} \begin{verbatim} (let ((vs '((1 1) (2 3) (4 5) (-1 2))) (2-norm (lizfcm.vector:p-norm 2))) @@ -124,7 +126,7 @@ x & y & 2norm\\[0pt] \end{center} \section{Question Four - |v|\textsubscript{1}} -\label{sec:org2196087} +\label{sec:org2b67b3e} \begin{verbatim} (let ((vs '((1 1) (2 3) (4 5) (-1 2))) (1-norm (lizfcm.vector:p-norm 1))) @@ -146,7 +148,7 @@ x & y & 1norm\\[0pt] \end{center} \section{Question Five - |v|\textsubscript{\(\infty\)}} -\label{sec:org11b8894} +\label{sec:org922206e} \begin{verbatim} (let ((vs '((1 1) (2 3) (4 5) (-1 2)))) (lizfcm.utils:table (:headers '("x" "y" "max-norm") @@ -167,11 +169,11 @@ x & y & infty-norm\\[0pt] \end{center} \section{Question Six - ||v - u|| via |v|\textsubscript{2}} -\label{sec:orga2324b2} +\label{sec:org29ec18f} \begin{verbatim} -(let* ((vs '((1 1) (2 3) (4 5) (-1 2))) - (vs2 '((7 9) (2 2) (8 -1) (4 4))) - (2-norm (lizfcm.vector:p-norm 2))) +(let ((vs '((1 1) (2 3) (4 5) (-1 2))) + (vs2 '((7 9) (2 2) (8 -1) (4 4))) + (2-norm (lizfcm.vector:p-norm 2))) (lizfcm.utils:table (:headers '("v1" "v2" "2-norm-d") :domain-order (v1 v2) :domain-values (mapcar (lambda (v1 v2) @@ -193,11 +195,11 @@ v1 & v2 & 2-norm\\[0pt] \end{center} \section{Question Seven - ||v - u|| via |v|\textsubscript{1}} -\label{sec:org388fbc7} +\label{sec:org7a87810} \begin{verbatim} -(let* ((vs '((1 1) (2 3) (4 5) (-1 2))) - (vs2 '((7 9) (2 2) (8 -1) (4 4))) - (1-norm (lizfcm.vector:p-norm 1))) +(let ((vs '((1 1) (2 3) (4 5) (-1 2))) + (vs2 '((7 9) (2 2) (8 -1) (4 4))) + (1-norm (lizfcm.vector:p-norm 1))) (lizfcm.utils:table (:headers '("v1" "v2" "1-norm-d") :domain-order (v1 v2) :domain-values (mapcar (lambda (v1 v2) @@ -219,10 +221,10 @@ v1 & v2 & 1-norm-d\\[0pt] \end{center} \section{Question Eight - ||v - u|| via |v|\textsubscript{\(\infty\)}} -\label{sec:org6e77f76} +\label{sec:org0f3b64f} \begin{verbatim} -(let* ((vs '((1 1) (2 3) (4 5) (-1 2))) - (vs2 '((7 9) (2 2) (8 -1) (4 4)))) +(let ((vs '((1 1) (2 3) (4 5) (-1 2))) + (vs2 '((7 9) (2 2) (8 -1) (4 4)))) (lizfcm.utils:table (:headers '("v1" "v2" "max-norm-d") :domain-order (v1 v2) :domain-values (mapcar (lambda (v1 v2) diff --git a/homeworks/hw-3.org b/homeworks/hw-3.org new file mode 100644 index 0000000..efeaa53 --- /dev/null +++ b/homeworks/hw-3.org @@ -0,0 +1,244 @@ +#+TITLE: HW 03 +#+AUTHOR: Elizabeth Hunt +#+STARTUP: entitiespretty fold inlineimages +#+LATEX_HEADER: \notindent \notag \usepackage{amsmath} \usepackage[a4paper,margin=1in,portrait]{geometry} +#+LATEX: \setlength\parindent{0pt} +#+OPTIONS: toc:nil + +* Question One +** Three Terms +\begin{align*} +Si_3(x) &= \int_0^x \frac{s - \frac{s^3}{3!} + \frac{s^5}{5!}}{s} dx \\ +&= x - \frac{x^3}{(3!)(3)} + \frac{x^5}{(5!)(5)} +\end{align*} +** Five Terms +\begin{align*} +Si_3(x) &= \int_0^x \frac{s - \frac{s^3}{3!} + \frac{s^5}{5!} - \frac{s^7}{7!} + \frac{s^9}{9!}}{s} dx \\ +&= x - \frac{x^3}{(3!)(3)} + \frac{x^5}{(5!)(5)} - \frac{x^7}{(7!)(7)} + \frac{s^9}{(9!)(9)} +\end{align*} +** Ten Terms +\begin{align*} +Si_{10}(x) &= \int_0^x \frac{s - \frac{s^3}{3!} + \frac{s^5}{5!} - \frac{s^7}{7!} + \frac{s^9}{9!} - \frac{s^{11}}{11!} + \frac{s^{13}}{13!} - \frac{s^{15}}{15!} + \frac{s^{17}}{17!} - \frac{s^{19}}{19!}}{s} ds \\ +&= x - \frac{x^3}{(3!)(3)} + \frac{x^5}{(5!)(5)} - \frac{x^7}{(7!)(7)} + \frac{s^9}{(9!)(9)} - \frac{s^{11}}{(11!)(11)} + \frac{s^{13}}{(13!)(13)} - \frac{s^{15}}{(15!)(15)} \\ +&+ \frac{s^{17}}{(17!)(17)} - \frac{s^{19}}{(19!)(19)} +\end{align*} +* Question Three +For the second term in the difference quotient, we can expand the taylor series centered at x=a: + +\begin{equation*} +f(x) = f(a) + f'(a)(x-a) + \frac{f''(a)}{2}(x-a)^2 + \cdots \\ +\end{equation*} + +Which we substitute into the difference quotient: + +\begin{equation*} +\frac{f(a) - f(a - h)}{h} = \frac{f(a) - (f(a) + f'(a)(x-a) + \frac{f''(a)}{2}(x-a)^2 + \cdots)}{h} +\end{equation*} + +And subs. $x=a-h$: + +\begin{align*} +\frac{f(a) - (f(a) + f'(a)(x-a) + \frac{f''(a)}{2}(x-a)^2 + \cdots)}{h} &= -f'(a)(-1) + -\frac{1}{2}f''(a)h \\ +&= f'(a) - \frac{1}{2}f''(a)h + \cdots \\ +\end{align*} + +Which we now plug into the initial $e_{\text{abs}}$: + +\begin{align*} +e_{\text{abs}} &= |f'(a) - \frac{f(a) - f(a - h)}{h}| \\ +&= |f'(a) - (f'(a) + -\frac{f''(a)}{2}h + \cdots)| \\ +&= |- \frac{1}{2}f''(a)h + \cdots | \\ +\end{align*} + +With the Taylor Remainder theorem we can absorb the series following the second term: + +\begin{equation*} +e_{\text{abs}} = |- \frac{1}{2}f''(a)h + \cdots | = |\frac{1}{2}f''(\xi)h| \leq Ch +\end{equation*} + +Thus our error is bounded linearly with $h$. + +* Question Four +For the first term in the difference quotient we know, from the given notes, + +\begin{equation*} +f(a+h) = f(a) + f'(a)h + \frac{1}{2}f''(a)h^2 + \frac{1}{6}f'''(a)(h^3) +\end{equation*} + +And from some of the work in Question Three, + +\begin{equation*} +f(a - h) = f(a) + f'(a)(-h) + \frac{1}{2}f''(a)(-h)^2 + \frac{1}{6}f'''(a)(-h^3) +\end{equation*} + +We can substitute immediately into $e_{\text{abs}} = |f'(a) - (\frac{f(a+h) - f(a-h)}{2h})|$: + +\begin{align*} +e_{\text{abs}} &= |f'(a) - \frac{1}{2h}((f(a) + f'(a)h + \frac{1}{2}f''(a)h^2 + \cdots) - (f(a) - f'(a)h + \frac{1}{2}f''(a)h^2 + \cdots))| \\ +&= |f'(a) - \frac{1}{2h}(2f'(a)h + \frac{1}{6}f'''(a)h^3 + \cdots)| \\ +&= |f'(a) - f'(a) - \frac{1}{12}f'''(a)h^2 + \cdots| \\ +&= |-\frac{1}{12}f'''(a)h^2 + \cdots| +\end{align*} + +Finally, with the Taylor Remainder theorem we can absorb the series following the third term: + +\begin{equation*} +e_{\text{abs}} = |-\frac{1}{12}f'''(\xi)h^2| = |\frac{1}{12}f'''(\xi)h^2| \leq Ch^2 +\end{equation*} + +Meaning that as $h$ scales linearly, our error is bounded by $h^2$ as opposed to linearly as in Question Three. + +* Question Six +** A +#+BEGIN_SRC lisp + (load "../lizfcm.asd") + (ql:quickload :lizfcm) + + (defun f (x) + (/ (- x 1) (+ x 1))) + + (defun fprime (x) + (/ 2 (expt (+ x 1) 2))) + + (let ((domain-values (loop for a from 0 to 2 + append + (loop for i from 0 to 9 + for h = (/ 1.0 (expt 2 i)) + collect (list a h))))) + (lizfcm.utils:table (:headers '("a" "h" "f'" "\\approx f'" "e_{\\text{abs}}") + :domain-order (a h) + :domain-values domain-values) + (fprime a) + (lizfcm.approx:fwd-derivative-at 'f a h) + (abs (- (fprime a) + (lizfcm.approx:fwd-derivative-at 'f a h))))) +#+END_SRC + +#+RESULTS: +| a | h | f' | \approx f' | e_{\text{abs}} | +| 0 | 1.0 | 2 | 1.0 | 1.0 | +| 0 | 0.5 | 2 | 1.3333333 | 0.66666675 | +| 0 | 0.25 | 2 | 1.5999999 | 0.4000001 | +| 0 | 0.125 | 2 | 1.7777777 | 0.22222233 | +| 0 | 0.0625 | 2 | 1.8823528 | 0.11764717 | +| 0 | 0.03125 | 2 | 1.939394 | 0.060606003 | +| 0 | 0.015625 | 2 | 1.9692307 | 0.030769348 | +| 0 | 0.0078125 | 2 | 1.9844971 | 0.01550293 | +| 0 | 0.00390625 | 2 | 1.992218 | 0.0077819824 | +| 0 | 0.001953125 | 2 | 1.9960938 | 0.00390625 | +| 1 | 1.0 | 1/2 | 0.33333334 | 0.16666666 | +| 1 | 0.5 | 1/2 | 0.4 | 0.099999994 | +| 1 | 0.25 | 1/2 | 0.44444445 | 0.055555552 | +| 1 | 0.125 | 1/2 | 0.47058824 | 0.029411763 | +| 1 | 0.0625 | 1/2 | 0.4848485 | 0.015151501 | +| 1 | 0.03125 | 1/2 | 0.4923077 | 0.0076923072 | +| 1 | 0.015625 | 1/2 | 0.49612403 | 0.0038759708 | +| 1 | 0.0078125 | 1/2 | 0.49805447 | 0.0019455254 | +| 1 | 0.00390625 | 1/2 | 0.49902534 | 0.00097465515 | +| 1 | 0.001953125 | 1/2 | 0.4995122 | 0.0004878044 | +| 2 | 1.0 | 2/9 | 0.16666666 | 0.055555567 | +| 2 | 0.5 | 2/9 | 0.19047618 | 0.031746045 | +| 2 | 0.25 | 2/9 | 0.2051282 | 0.017094031 | +| 2 | 0.125 | 2/9 | 0.21333337 | 0.008888856 | +| 2 | 0.0625 | 2/9 | 0.21768713 | 0.004535094 | +| 2 | 0.03125 | 2/9 | 0.21993065 | 0.002291575 | +| 2 | 0.015625 | 2/9 | 0.22106934 | 0.0011528879 | +| 2 | 0.0078125 | 2/9 | 0.22164536 | 0.00057686865 | +| 2 | 0.00390625 | 2/9 | 0.22193146 | 0.00029076636 | +| 2 | 0.001953125 | 2/9 | 0.22207642 | 0.00014580786 | + + +* Question Nine +** C + +#+BEGIN_SRC lisp + (load "../lizfcm.asd") + (ql:quickload :lizfcm) + + (defun factorial (n) + (if (= n 0) + 1 + (* n (factorial (- n 1))))) + + (defun taylor-term (n x) + (/ (* (expt (- 1) n) + (expt x (+ (* 2 n) 1))) + (* (factorial n) + (+ (* 2 n) 1)))) + + (defun f (x &optional (max-iterations 30)) + (let ((sum 0.0)) + (dotimes (n max-iterations) + (setq sum (+ sum (taylor-term n x)))) + (* sum (/ 2 (sqrt pi))))) + + (defun fprime (x) + (* (/ 2 (sqrt pi)) (exp (- 0 (* x x))))) + + (let ((domain-values (loop for a from 0 to 1 + append + (loop for i from 0 to 9 + for h = (/ 1.0 (expt 2 i)) + collect (list a h))))) + (lizfcm.utils:table (:headers '("a" "h" "f'" "\\approx f'" "e_{\\text{abs}}") + :domain-order (a h) + :domain-values domain-values) + (fprime a) + (lizfcm.approx:central-derivative-at 'f a h) + (abs (- (fprime a) + (lizfcm.approx:central-derivative-at 'f a h))))) +#+END_SRC + + +| a | h | f' | \approx f' | e_{\text{abs}} | +| 0 | 1.0 | 1.1283791670955126d0 | 0.8427006725464232d0 | 0.28567849454908933d0 | +| 0 | 0.5 | 1.1283791670955126d0 | 1.0409997446922075d0 | 0.0873794224033051d0 | +| 0 | 0.25 | 1.1283791670955126d0 | 1.1053055663206806d0 | 0.023073600774832004d0 | +| 0 | 0.125 | 1.1283791670955126d0 | 1.122529655394656d0 | 0.005849511700856569d0 | +| 0 | 0.0625 | 1.1283791670955126d0 | 1.1269116944798618d0 | 0.0014674726156507223d0 | +| 0 | 0.03125 | 1.1283791670955126d0 | 1.1280120131008824d0 | 3.6715399463016496d-4 | +| 0 | 0.015625 | 1.1283791670955126d0 | 1.1282873617826952d0 | 9.180531281738347d-5 | +| 0 | 0.0078125 | 1.1283791670955126d0 | 1.128356232581468d0 | 2.293451404455915d-5 | +| 0 | 0.00390625 | 1.1283791670955126d0 | 1.1283734502811613d0 | 5.71681435124205d-6 | +| 0 | 0.001953125 | 1.1283791670955126d0 | 1.1283777547060847d0 | 1.4123894278572635d-6 | +| 1 | 1.0 | 0.41510750774498784d0 | 0.4976611317561498d0 | 0.08255362401116195d0 | +| 1 | 0.5 | 0.41510750774498784d0 | 0.44560523266293384d0 | 0.030497724917946d0 | +| 1 | 0.25 | 0.41510750774498784d0 | 0.4234889628937013d0 | 0.008381455148713468d0 | +| 1 | 0.125 | 0.41510750774498784d0 | 0.41725265825950153d0 | 0.002145150514513694d0 | +| 1 | 0.0625 | 0.41510750774498784d0 | 0.41564710776310854d0 | 5.396000181207006d-4 | +| 1 | 0.03125 | 0.41510750774498784d0 | 0.4152414157140871d0 | 1.3390796909928948d-4 | +| 1 | 0.015625 | 0.41510750774498784d0 | 0.41514241394084905d0 | 3.490619586121735d-5 | +| 1 | 0.0078125 | 0.41510750774498784d0 | 0.41510582632900395d0 | 1.6814159838896003d-6 | +| 1 | 0.00390625 | 0.41510750774498784d0 | 0.415092913054238d0 | 1.4594690749825112d-5 | +| 1 | 0.001953125 | 0.41510750774498784d0 | 0.4150670865046777d0 | 4.0421240310117845d-5 | + + +* Question Twelve + +First we'll place a bound on $h$; looking at a graph of $f$ it's pretty obvious from the asymptotes that we don't want to go much further than $|h| = 2 - \frac{pi}{2}$. + +Following similar reasoning as Question Four, we can determine an optimal $h$ by computing $e_{\text{abs}}$ for the central difference, but now including a roundoff error for each time we run $f$ +such that $|f_{\text{machine}}(x) - f(x)| \le \epsilon_{\text{dblprec}}$ (we'll use double precision numbers, from HW 2 we know $\epsilon_{\text{dblprec}} \approx 2.22045 (10^{-16})$). + +We'll just assume $|f_{\text{machine}}(x) - f(x)| = \epsilon_{\text{dblprec}}$ so our new difference quotient becomes: + +\begin{align*} +e_{\text{abs}} &= |f'(a) - (\frac{f(a+h) - f(a-h) + 2\epsilon_{\text{dblprec}}}{2h})| \\ +&= |\frac{1}{12}f'''(\xi)h^2 + \frac{\epsilon_{\text{dblprec}}}{h}| +\end{align*} + +Because we bounded our $|h| = 2 - \frac{pi}{2}$ we'll find the maximum value of $f'''$ between $a - (2 - \frac{\pi}{2})$ and $a - (2 - \frac{\pi}{3})$. Using [[https://www.desmos.com/calculator/gen1zpohh2][desmos]] I found this to be -2. + +Thus, $e_{\text{abs}} \leq \frac{1}{6}h^2 + \frac{\epsilon_{\text{dblprec}}}{h}$. Finding the derivative: + +\begin{equation*} +e' = \frac{1}{3}h - \frac{\epsilon_{\text{dblprec}}}{h^2} +\end{equation*} + +And solving at $e' = 0$: + +\begin{equation*} +\frac{1}{3}h = \frac{\epsilon_{\text{dblprec}}}{h^2} \Rightarrow h^3 = 3\epsilon_{\text{dblprec}} \Rightarrow h = (3\epsilon_{\text{dblprec}})^{1/3} +\end{equation*} + +Which is $\approx (3(2.22045 (10^{-16}))^{\frac{1}{3}} \approx 8.7335 10^{-6}$. diff --git a/homeworks/hw-3.pdf b/homeworks/hw-3.pdf Binary files differnew file mode 100644 index 0000000..1f9bac6 --- /dev/null +++ b/homeworks/hw-3.pdf diff --git a/homeworks/hw-3.tex b/homeworks/hw-3.tex new file mode 100644 index 0000000..b3d029d --- /dev/null +++ b/homeworks/hw-3.tex @@ -0,0 +1,250 @@ +% Created 2023-10-07 Sat 14:49 +% Intended LaTeX compiler: pdflatex +\documentclass[11pt]{article} +\usepackage[utf8]{inputenc} +\usepackage[T1]{fontenc} +\usepackage{graphicx} +\usepackage{longtable} +\usepackage{wrapfig} +\usepackage{rotating} +\usepackage[normalem]{ulem} +\usepackage{amsmath} +\usepackage{amssymb} +\usepackage{capt-of} +\usepackage{hyperref} +\notindent \notag \usepackage{amsmath} \usepackage[a4paper,margin=1in,portrait]{geometry} +\author{Elizabeth Hunt} +\date{\today} +\title{HW 03} +\hypersetup{ + pdfauthor={Elizabeth Hunt}, + pdftitle={HW 03}, + pdfkeywords={}, + pdfsubject={}, + pdfcreator={Emacs 28.2 (Org mode 9.7-pre)}, + pdflang={English}} +\begin{document} + +\maketitle +\setlength\parindent{0pt} + +\section{Question One} +\label{sec:org6f2bd27} +\subsection{Three Terms} +\label{sec:orgeb827ff} +\begin{align*} +Si_3(x) &= \int_0^x \frac{s - \frac{s^3}{3!} + \frac{s^5}{5!}}{s} dx \\ +&= x - \frac{x^3}{(3!)(3)} + \frac{x^5}{(5!)(5)} +\end{align*} +\subsection{Five Terms} +\label{sec:orge6a15e4} +\begin{align*} +Si_3(x) &= \int_0^x \frac{s - \frac{s^3}{3!} + \frac{s^5}{5!} - \frac{s^7}{7!} + \frac{s^9}{9!}}{s} dx \\ +&= x - \frac{x^3}{(3!)(3)} + \frac{x^5}{(5!)(5)} - \frac{x^7}{(7!)(7)} + \frac{s^9}{(9!)(9)} +\end{align*} +\subsection{Ten Terms} +\label{sec:orge87e346} +\begin{align*} +Si_{10}(x) &= \int_0^x \frac{s - \frac{s^3}{3!} + \frac{s^5}{5!} - \frac{s^7}{7!} + \frac{s^9}{9!} - \frac{s^{11}}{11!} + \frac{s^{13}}{13!} - \frac{s^{15}}{15!} + \frac{s^{17}}{17!} - \frac{s^{19}}{19!}}{s} ds \\ +&= x - \frac{x^3}{(3!)(3)} + \frac{x^5}{(5!)(5)} - \frac{x^7}{(7!)(7)} + \frac{s^9}{(9!)(9)} - \frac{s^{11}}{(11!)(11)} + \frac{s^{13}}{(13!)(13)} - \frac{s^{15}}{(15!)(15)} \\ +&+ \frac{s^{17}}{(17!)(17)} - \frac{s^{19}}{(19!)(19)} +\end{align*} +\section{Question Three} +\label{sec:org6e2f7fc} +For the second term in the difference quotient, we can expand the taylor series centered at x=a: + +\begin{equation*} +f(x) = f(a) + f'(a)(x-a) + \frac{f''(a)}{2}(x-a)^2 + \cdots \\ +\end{equation*} + +Which we substitute into the difference quotient: + +\begin{equation*} +\frac{f(a) - f(a - h)}{h} = \frac{f(a) - (f(a) + f'(a)(x-a) + \frac{f''(a)}{2}(x-a)^2 + \cdots)}{h} +\end{equation*} + +And subs. \(x=a-h\): + +\begin{align*} +\frac{f(a) - (f(a) + f'(a)(x-a) + \frac{f''(a)}{2}(x-a)^2 + \cdots)}{h} &= -f'(a)(-1) + -\frac{1}{2}f''(a)h \\ +&= f'(a) - \frac{1}{2}f''(a)h + \cdots \\ +\end{align*} + +Which we now plug into the initial \(e_{\text{abs}}\): + +\begin{align*} +e_{\text{abs}} &= |f'(a) - \frac{f(a) - f(a - h)}{h}| \\ +&= |f'(a) - (f'(a) + -\frac{f''(a)}{2}h + \cdots)| \\ +&= |- \frac{1}{2}f''(a)h + \cdots | \\ +\end{align*} + +With the Taylor Remainder theorem we can absorb the series following the second term: + +\begin{equation*} +e_{\text{abs}} = |- \frac{1}{2}f''(a)h + \cdots | = |\frac{1}{2}f''(\xi)h| \leq Ch +\end{equation*} + +Thus our error is bounded linearly with \(h\). + +\section{Question Four} +\label{sec:orga7d02a2} +For the first term in the difference quotient we know, from the given notes, + +\begin{equation*} +f(a+h) = f(a) + f'(a)h + \frac{1}{2}f''(a)h^2 + \frac{1}{6}f'''(a)(h^3) +\end{equation*} + +And from some of the work in Question Three, + +\begin{equation*} +f(a - h) = f(a) + f'(a)(-h) + \frac{1}{2}f''(a)(-h)^2 + \frac{1}{6}f'''(a)(-h^3) +\end{equation*} + +We can substitute immediately into \(e_{\text{abs}} = |f'(a) - (\frac{f(a+h) - f(a-h)}{2h})|\): + +\begin{align*} +e_{\text{abs}} &= |f'(a) - \frac{1}{2h}((f(a) + f'(a)h + \frac{1}{2}f''(a)h^2 + \cdots) - (f(a) - f'(a)h + \frac{1}{2}f''(a)h^2 + \cdots))| \\ +&= |f'(a) - \frac{1}{2h}(2f'(a)h + \frac{1}{6}f'''(a)h^3 + \cdots)| \\ +&= |f'(a) - f'(a) - \frac{1}{12}f'''(a)h^2 + \cdots| \\ +&= |-\frac{1}{12}f'''(a)h^2 + \cdots| +\end{align*} + +Finally, with the Taylor Remainder theorem we can absorb the series following the third term: + +\begin{equation*} +e_{\text{abs}} = |-\frac{1}{12}f'''(\xi)h^2| = |\frac{1}{12}f'''(\xi)h^2| \leq Ch^2 +\end{equation*} + +Meaning that as \(h\) scales linearly, our error is bounded by \(h^2\) as opposed to linearly as in Question Three. + +\section{Question Six} +\label{sec:org7b05811} +\subsection{A} +\label{sec:org8341a77} +\begin{verbatim} +(load "../lizfcm.asd") +(ql:quickload :lizfcm) + +(defun f (x) + (/ (- x 1) (+ x 1))) + +(defun fprime (x) + (/ 2 (expt (+ x 1) 2))) + +(let ((domain-values (loop for a from 0 to 2 + append + (loop for i from 0 to 9 + for h = (/ 1.0 (expt 2 i)) + collect (list a h))))) + (lizfcm.utils:table (:headers '("a" "h" "f'" "\\approx f'" "e_{\\text{abs}}") + :domain-order (a h) + :domain-values domain-values) + (fprime a) + (lizfcm.approx:fwd-derivative-at 'f a h) + (abs (- (fprime a) + (lizfcm.approx:fwd-derivative-at 'f a h))))) +\end{verbatim} + + +\section{Question Nine} +\label{sec:orgeb1839f} +\subsection{C} +\label{sec:org5691277} + +\begin{verbatim} +(load "../lizfcm.asd") +(ql:quickload :lizfcm) + +(defun factorial (n) + (if (= n 0) + 1 + (* n (factorial (- n 1))))) + +(defun taylor-term (n x) + (/ (* (expt (- 1) n) + (expt x (+ (* 2 n) 1))) + (* (factorial n) + (+ (* 2 n) 1)))) + +(defun f (x &optional (max-iterations 30)) + (let ((sum 0.0)) + (dotimes (n max-iterations) + (setq sum (+ sum (taylor-term n x)))) + (* sum (/ 2 (sqrt pi))))) + +(defun fprime (x) + (* (/ 2 (sqrt pi)) (exp (- 0 (* x x))))) + +(let ((domain-values (loop for a from 0 to 1 + append + (loop for i from 0 to 9 + for h = (/ 1.0 (expt 2 i)) + collect (list a h))))) + (lizfcm.utils:table (:headers '("a" "h" "f'" "\\approx f'" "e_{\\text{abs}}") + :domain-order (a h) + :domain-values domain-values) + (fprime a) + (lizfcm.approx:central-derivative-at 'f a h) + (abs (- (fprime a) + (lizfcm.approx:central-derivative-at 'f a h))))) +\end{verbatim} + + +\begin{center} +\begin{tabular}{rrrrr} +a & h & f' & \(\approx\) f' & e\textsubscript{\text{abs}}\\[0pt] +0 & 1.0 & 1.1283791670955126d0 & 0.8427006725464232d0 & 0.28567849454908933d0\\[0pt] +0 & 0.5 & 1.1283791670955126d0 & 1.0409997446922075d0 & 0.0873794224033051d0\\[0pt] +0 & 0.25 & 1.1283791670955126d0 & 1.1053055663206806d0 & 0.023073600774832004d0\\[0pt] +0 & 0.125 & 1.1283791670955126d0 & 1.122529655394656d0 & 0.005849511700856569d0\\[0pt] +0 & 0.0625 & 1.1283791670955126d0 & 1.1269116944798618d0 & 0.0014674726156507223d0\\[0pt] +0 & 0.03125 & 1.1283791670955126d0 & 1.1280120131008824d0 & 3.6715399463016496d-4\\[0pt] +0 & 0.015625 & 1.1283791670955126d0 & 1.1282873617826952d0 & 9.180531281738347d-5\\[0pt] +0 & 0.0078125 & 1.1283791670955126d0 & 1.128356232581468d0 & 2.293451404455915d-5\\[0pt] +0 & 0.00390625 & 1.1283791670955126d0 & 1.1283734502811613d0 & 5.71681435124205d-6\\[0pt] +0 & 0.001953125 & 1.1283791670955126d0 & 1.1283777547060847d0 & 1.4123894278572635d-6\\[0pt] +1 & 1.0 & 0.41510750774498784d0 & 0.4976611317561498d0 & 0.08255362401116195d0\\[0pt] +1 & 0.5 & 0.41510750774498784d0 & 0.44560523266293384d0 & 0.030497724917946d0\\[0pt] +1 & 0.25 & 0.41510750774498784d0 & 0.4234889628937013d0 & 0.008381455148713468d0\\[0pt] +1 & 0.125 & 0.41510750774498784d0 & 0.41725265825950153d0 & 0.002145150514513694d0\\[0pt] +1 & 0.0625 & 0.41510750774498784d0 & 0.41564710776310854d0 & 5.396000181207006d-4\\[0pt] +1 & 0.03125 & 0.41510750774498784d0 & 0.4152414157140871d0 & 1.3390796909928948d-4\\[0pt] +1 & 0.015625 & 0.41510750774498784d0 & 0.41514241394084905d0 & 3.490619586121735d-5\\[0pt] +1 & 0.0078125 & 0.41510750774498784d0 & 0.41510582632900395d0 & 1.6814159838896003d-6\\[0pt] +1 & 0.00390625 & 0.41510750774498784d0 & 0.415092913054238d0 & 1.4594690749825112d-5\\[0pt] +1 & 0.001953125 & 0.41510750774498784d0 & 0.4150670865046777d0 & 4.0421240310117845d-5\\[0pt] +\end{tabular} +\end{center} + + +\section{Question Twelve} +\label{sec:orgc55bfd1} + +First we'll place a bound on \(h\); looking at a graph of \(f\) it's pretty obvious from the asymptotes that we don't want to go much further than \(|h| = 2 - \frac{pi}{2}\). + +Following similar reasoning as Question Four, we can determine an optimal \(h\) by computing \(e_{\text{abs}}\) for the central difference, but now including a roundoff error for each time we run \(f\) +such that \(|f_{\text{machine}}(x) - f(x)| \le \epsilon_{\text{dblprec}}\) (we'll use double precision numbers, from HW 2 we know \(\epsilon_{\text{dblprec}} \approx 2.22045 (10^{-16})\)). + +We'll just assume \(|f_{\text{machine}}(x) - f(x)| = \epsilon_{\text{dblprec}}\) so our new difference quotient becomes: + +\begin{align*} +e_{\text{abs}} &= |f'(a) - (\frac{f(a+h) - f(a-h) + 2\epsilon_{\text{dblprec}}}{2h})| \\ +&= |\frac{1}{12}f'''(\xi)h^2 + \frac{\epsilon_{\text{dblprec}}}{h}| +\end{align*} + +Because we bounded our \(|h| = 2 - \frac{pi}{2}\) we'll find the maximum value of \(f'''\) between \(a - (2 - \frac{\pi}{2})\) and \(a - (2 - \frac{\pi}{3})\). Using \href{https://www.desmos.com/calculator/gen1zpohh2}{desmos} I found this to be -2. + +Thus, \(e_{\text{abs}} \leq \frac{1}{6}h^2 + \frac{\epsilon_{\text{dblprec}}}{h}\). Finding the derivative: + +\begin{equation*} +e' = \frac{1}{3}h - \frac{\epsilon_{\text{dblprec}}}{h^2} +\end{equation*} + +And solving at \(e' = 0\): + +\begin{equation*} +\frac{1}{3}h = \frac{\epsilon_{\text{dblprec}}}{h^2} \Rightarrow h^3 = 3\epsilon_{\text{dblprec}} \Rightarrow h = (3\epsilon_{\text{dblprec}})^{1/3} +\end{equation*} + +Which is \(\approx (3(2.22045 (10^{-16}))^{\frac{1}{3}} \approx 8.7335 10^{-6}\). +\end{document}
\ No newline at end of file diff --git a/notes/Oct-4.org b/notes/Oct-4.org new file mode 100644 index 0000000..8b8466f --- /dev/null +++ b/notes/Oct-4.org @@ -0,0 +1,22 @@ +[[ a_{11} a_{12} \cdots a_{1n} | b_1] + [ 0 (a_{22} - \frac{a_{}_{21}}{a_{22}}a_{11}) \cdots a_{2n} | b_2 - \frac{a_{21}}{a_{11}}b_1 ]] + +#+BEGIN_SRC c + for (int i = 1; i < n; i++) { + float factor = -a[i][0] / a[0][0]; + for (int j = 1; j < n; j++) { + a[i][j] = a[i][j] + factor * a[0][j]; + } + b[i] = b[i] + factor * b[0]; + } + + for (int k = 0; k < (n - 1); k++) { + for (int i = k+1; i < n; i++) { + float factor = -a[i][k] / a[k][k]; + for (int j = k+1; j < n; j++) { + a[i][j] = a[i][j] + factor * a[j][k]; + } + b[i] = b[i] + factor*b[k]; + } + } +#+END_SRC diff --git a/notes/Oct-6.org b/notes/Oct-6.org new file mode 100644 index 0000000..8cbff29 --- /dev/null +++ b/notes/Oct-6.org @@ -0,0 +1,13 @@ +#+BEGIN_SRC c + for (int k = 0; i < (n - 1); k++) { + for (int i = k+1; i< n; i++) { + float factor = a[i][k] / a[k][k]; + for (int j = k+1; j < k; j++) { + a[i][j] = a[i][j] - factor * a[k][j]; + } + b[i] = b[i] - factor * b[k]; + } + } +#+END_SRC + + diff --git a/notes/Sep-15.org b/notes/Sep-15.org index d5bf371..8a64089 100644 --- a/notes/Sep-15.org +++ b/notes/Sep-15.org @@ -36,11 +36,11 @@ Again, $f'(a) \approx \frac{f(a+h) - f(a)}{h}$, $e = |\frac{1}{2} f''(a) + \frac{1}{3!}h^2 f'''(a) + \cdots$ -$R_2 = \frac{h}{2} f''(u)$ +$R_2 = \frac{h}{2} f''(\xi)$ -$|\frac{h}{2} f''(u)| \leq M h^1$ +$|\frac{h}{2} f''(\xi)| \leq M h^1$ -$M = \frac{1}{2}|f'(u)|$ +$M = \frac{1}{2}|f'(\xi)|$ *** Another approximation @@ -48,5 +48,5 @@ $\text{err} = |f'(a) - \frac{f(a) - f(a - h)}{h}|$ $= f'(a) - \frac{1}{h}(f(a) - (f(a) + f'(a)(a - (a - h)) + \frac{1}{2}f''(a)(a-(a-h))^2 + \cdots))$ -$= |f'(a) - \frac{1}{h}(f'(a) + \frac{1}{2}f''(a)h)|$ +$= |f'(a) - (f'(a) + \frac{1}{2}f''(a)h)|$ diff --git a/notes/Sep-15.pdf b/notes/Sep-15.pdf Binary files differnew file mode 100644 index 0000000..43a4f34 --- /dev/null +++ b/notes/Sep-15.pdf diff --git a/notes/Sep-15.tex b/notes/Sep-15.tex new file mode 100644 index 0000000..52610ed --- /dev/null +++ b/notes/Sep-15.tex @@ -0,0 +1,88 @@ +% Created 2023-09-29 Fri 10:00 +% Intended LaTeX compiler: pdflatex +\documentclass[11pt]{article} +\usepackage[utf8]{inputenc} +\usepackage[T1]{fontenc} +\usepackage{graphicx} +\usepackage{longtable} +\usepackage{wrapfig} +\usepackage{rotating} +\usepackage[normalem]{ulem} +\usepackage{amsmath} +\usepackage{amssymb} +\usepackage{capt-of} +\usepackage{hyperref} +\author{Elizabeth Hunt} +\date{\today} +\title{} +\hypersetup{ + pdfauthor={Elizabeth Hunt}, + pdftitle={}, + pdfkeywords={}, + pdfsubject={}, + pdfcreator={Emacs 28.2 (Org mode 9.7-pre)}, + pdflang={English}} +\begin{document} + +\tableofcontents + +\section{Taylor Series Approx.} +\label{sec:orgcc72ed1} +Suppose f has \(\infty\) many derivatives near a point a. Then the taylor series is given by + +\(f(x) = \Sigma_{n=0}^{\infty} \frac{f^{(n)}(a)}{n!}(x-a)^n\) + +For increment notation we can write + +\(f(a + h) = f(a) + f'(a)(a+h - a) + \dots\) + +\(= \Sigma_{n=0}^{\infty} \frac{f^{(n)}(a)}{h!} (h^n)\) + +Consider the approximation + +\(e = |f'(a) - \frac{f(a + h) - f(a)}{h}| = |f'(a) - \frac{1}{h}(f(a + h) - f(a))|\) + +Substituting\ldots{} + +\(= |f'(a) - \frac{1}{h}((f(a) + f'(a) h + \frac{f''(a)}{2} h^2 + \cdots) - f(a))|\) + +\(f(a) - f(a) = 0\)\ldots{} and \(distribute the h\) + +\(= |-1/2 f''(a) h + \frac{1}{6}f'''(a)h^2 \cdots|\) + +\subsection{With Remainder} +\label{sec:org7dfd6c7} +We can determine for some u \(f(a + h) = f(a) + f'(a)h + \frac{1}{2}f''(u)h^2\) + +and so the error is \(e = |f'(a) - \frac{f(a + h) - f(a)}{h}| = |\frac{h}{2}f''(u)|\) + +\begin{itemize} +\item\relax [\url{https://openstax.org/books/calculus-volume-2/pages/6-3-taylor-and-maclaurin-series}] +\begin{itemize} +\item > Taylor's Theorem w/ Remainder +\end{itemize} +\end{itemize} + + +\subsection{Of Deriviatives} +\label{sec:org1ec7c9b} + +Again, \(f'(a) \approx \frac{f(a+h) - f(a)}{h}\), + +\(e = |\frac{1}{2} f''(a) + \frac{1}{3!}h^2 f'''(a) + \cdots\) + +\(R_2 = \frac{h}{2} f''(u)\) + +\(|\frac{h}{2} f''(u)| \leq M h^1\) + +\(M = \frac{1}{2}|f'(u)|\) + +\subsubsection{Another approximation} +\label{sec:org16193b9} + +\(\text{err} = |f'(a) - \frac{f(a) - f(a - h)}{h}|\) + +\(= f'(a) - \frac{1}{h}(f(a) - (f(a) + f'(a)(a - (a - h)) + \frac{1}{2}f''(a)(a-(a-h))^2 + \cdots))\) + +\(= |f'(a) - \frac{1}{h}(f'(a) + \frac{1}{2}f''(a)h)|\) +\end{document}
\ No newline at end of file diff --git a/src/approx/derivative.lisp b/src/approx/derivative.lisp index 02fcb4c..8aa171a 100644 --- a/src/approx/derivative.lisp +++ b/src/approx/derivative.lisp @@ -1,9 +1,17 @@ (in-package :lizfcm.approx) -(defun derivative-at (f x &optional (delta 0.01)) +(defun central-derivative-at (f x &optional (delta 0.01)) (let* ((x2 (+ x delta)) (x1 (- x delta)) (y2 (apply f (list x2))) (y1 (apply f (list x1)))) (/ (- y2 y1) (- x2 x1)))) + +(defun fwd-derivative-at (f x &optional (delta 0.01)) + (let* ((x2 (+ x delta)) + (x1 x) + (y2 (apply f (list x2))) + (y1 (apply f (list x1)))) + (/ (- y2 y1) + (- x2 x1)))) diff --git a/src/approx/maceps.lisp b/src/approx/maceps.lisp index debaa67..e2738e4 100644 --- a/src/approx/maceps.lisp +++ b/src/approx/maceps.lisp @@ -3,10 +3,10 @@ (defun compute-maceps (f a init) (let ((h init) (err init)) - (loop while (> err 0) + (loop collect (list a h err) do (setf h (/ h 2) err (abs (- (funcall f (+ a h)) (funcall f a)))) - when (> err 0) - collect (list a h err)))) + while (> err 0)))) + diff --git a/src/approx/package.lisp b/src/approx/package.lisp index c11eb60..67c6a90 100644 --- a/src/approx/package.lisp +++ b/src/approx/package.lisp @@ -1,5 +1,6 @@ (in-package :cl-user) (defpackage lizfcm.approx (:use :cl) - (:export :derivative-at + (:export :central-derivative-at + :fwd-derivative-at :compute-maceps)) diff --git a/tests/approx.lisp b/tests/approx.lisp index 2fd8124..588b16d 100644 --- a/tests/approx.lisp +++ b/tests/approx.lisp @@ -11,7 +11,7 @@ :in lizfcm-test-suite) (in-suite approx-suite) -(test derivative-at +(test central-derivative-at :description "derivative at is within bounds" (let ((f (lambda (x) (* x x))) (x 2) @@ -19,6 +19,18 @@ (f-prime-at-x 4) (delta 0.01)) (is (within-range-p - (derivative-at f x delta) + (central-derivative-at f x delta) + f-prime-at-x + accepted-delta)))) + +(test fwd-derivative-at + :description "forward derivative at is within bounds" + (let ((f (lambda (x) (* x x))) + (x 2) + (accepted-delta 0.02) + (f-prime-at-x 4) + (delta 0.01)) + (is (within-range-p + (forward-derivative-at f x delta) f-prime-at-x accepted-delta)))) |