A. Lambert W function
In A1, you were asked to prove that the Lambert W function (also called the product-log function ) is a well-defined function W : [0,\infty) \to [0,\infty) satisfying W(a) e^{W(a)} = a . Define \Omega := W(1) . Last time, essentially by the change of sign theorem, you showed that \Omega \in [0.56, 0.58] .
β Show that \Omega is a fixed point of g and h as defined by
\begin{align}
g(x) = e^{-x}, \qquad h(x) = \log \left|\frac{1}{x}\right|
\end{align}
Answer.
We have \Omega = W(1) and so \Omega e^{\Omega} = 1 and \Omega = e^{-\Omega} = g(\Omega) . We also have e^{\Omega} = \Omega^{-1} = |\Omega^{-1}| because \Omega > 0 . Taking logarithms on both sides gives \Omega = \log |\Omega^{-1}| = h(\Omega) .
π» Show numerically that x_{n+1} = g(x_n) converges and y_{n+1} = h(y_n) diverges for a few different choices of x_1 and y_1 .
Answer.
We will use the function simple_iteration( g, x1 ) as defined above:
function printErrors ( x ) = [1 : length (x) x (@.abs (x- Ξ©))]# you can also add colour! = Highlighter (-> (j == 3 ) && (data[i, j] > 1e-2 ),"red bold" = Highlighter (-> (j == 3 ) && (data[i, j] < 1e-5 ),"green bold" pretty_table ( = ["n" , "x[n]" , "|x[n] - Ξ©|" ],= (h1, h2) end println ("Iteration xβββ = g(xβ) with xβ = 1" )printErrors ( simple_iteration ( g, 1 . ) )
Iteration xβββ = g(xβ) with xβ = 1
ββββββββ¬βββββββββββ¬ββββββββββββββ
β n β x[n] β |x[n] - Ξ©| β
ββββββββΌβββββββββββΌββββββββββββββ€
β 1.0 β 1.0 β 0.432857 β
β 2.0 β 0.367879 β 0.199264 β
β 3.0 β 0.692201 β 0.125057 β
β 4.0 β 0.500474 β 0.0666698 β
β 5.0 β 0.606244 β 0.0391002 β
β 6.0 β 0.545396 β 0.0217475 β
β 7.0 β 0.579612 β 0.012469 β
β 8.0 β 0.560115 β 0.00702783 β
β 9.0 β 0.571143 β 0.00399982 β
β 10.0 β 0.564879 β 0.00226394 β
β 11.0 β 0.568429 β 0.00128543 β
β 12.0 β 0.566415 β 0.000728557 β
β 13.0 β 0.567557 β 0.000413347 β
β 14.0 β 0.566909 β 0.000234378 β
β 15.0 β 0.567276 β 0.000132942 β
β 16.0 β 0.567068 β 7.5392e-5 β
β 17.0 β 0.567186 β 4.27597e-5 β
β 18.0 β 0.567119 β 2.42504e-5 β
β 19.0 β 0.567157 β 1.37536e-5 β
β 20.0 β 0.567135 β 7.8002e-6 β
β 21.0 β 0.567148 β 4.42385e-6 β
β 22.0 β 0.567141 β 2.50895e-6 β
β 23.0 β 0.567145 β 1.42294e-6 β
β 24.0 β 0.567142 β 8.07008e-7 β
β 25.0 β 0.567144 β 4.5769e-7 β
β 26.0 β 0.567143 β 2.59576e-7 β
β 27.0 β 0.567143 β 1.47217e-7 β
β 28.0 β 0.567143 β 8.34929e-8 β
β 29.0 β 0.567143 β 4.73524e-8 β
β 30.0 β 0.567143 β 2.68556e-8 β
β 31.0 β 0.567143 β 1.5231e-8 β
β 32.0 β 0.567143 β 8.63815e-9 β
β 33.0 β 0.567143 β 4.89907e-9 β
β 34.0 β 0.567143 β 2.77847e-9 β
β 35.0 β 0.567143 β 1.57579e-9 β
β 36.0 β 0.567143 β 8.93701e-10 β
β 37.0 β 0.567143 β 5.06855e-10 β
β 38.0 β 0.567143 β 2.87461e-10 β
β 39.0 β 0.567143 β 1.6303e-10 β
β 40.0 β 0.567143 β 9.24626e-11 β
β 41.0 β 0.567143 β 5.24385e-11 β
ββββββββ΄βββββββββββ΄ββββββββββββββ
= 100 .;println ("Iteration xβββ = g(xβ) with xβ = $ x1" )printErrors ( simple_iteration ( g, x1 ) )
Iteration xβββ = g(xβ) with xβ = 100.0
ββββββββ¬ββββββββββββββ¬ββββββββββββββ
β n β x[n] β |x[n] - Ξ©| β
ββββββββΌββββββββββββββΌββββββββββββββ€
β 1.0 β 100.0 β 99.4329 β
β 2.0 β 3.72008e-44 β 0.567143 β
β 3.0 β 1.0 β 0.432857 β
β 4.0 β 0.367879 β 0.199264 β
β 5.0 β 0.692201 β 0.125057 β
β 6.0 β 0.500474 β 0.0666698 β
β 7.0 β 0.606244 β 0.0391002 β
β 8.0 β 0.545396 β 0.0217475 β
β 9.0 β 0.579612 β 0.012469 β
β 10.0 β 0.560115 β 0.00702783 β
β 11.0 β 0.571143 β 0.00399982 β
β 12.0 β 0.564879 β 0.00226394 β
β 13.0 β 0.568429 β 0.00128543 β
β 14.0 β 0.566415 β 0.000728557 β
β 15.0 β 0.567557 β 0.000413347 β
β 16.0 β 0.566909 β 0.000234378 β
β 17.0 β 0.567276 β 0.000132942 β
β 18.0 β 0.567068 β 7.5392e-5 β
β 19.0 β 0.567186 β 4.27597e-5 β
β 20.0 β 0.567119 β 2.42504e-5 β
β 21.0 β 0.567157 β 1.37536e-5 β
β 22.0 β 0.567135 β 7.8002e-6 β
β 23.0 β 0.567148 β 4.42385e-6 β
β 24.0 β 0.567141 β 2.50895e-6 β
β 25.0 β 0.567145 β 1.42294e-6 β
β 26.0 β 0.567142 β 8.07008e-7 β
β 27.0 β 0.567144 β 4.5769e-7 β
β 28.0 β 0.567143 β 2.59576e-7 β
β 29.0 β 0.567143 β 1.47217e-7 β
β 30.0 β 0.567143 β 8.34929e-8 β
β 31.0 β 0.567143 β 4.73524e-8 β
β 32.0 β 0.567143 β 2.68556e-8 β
β 33.0 β 0.567143 β 1.5231e-8 β
β 34.0 β 0.567143 β 8.63815e-9 β
β 35.0 β 0.567143 β 4.89907e-9 β
β 36.0 β 0.567143 β 2.77847e-9 β
β 37.0 β 0.567143 β 1.57579e-9 β
β 38.0 β 0.567143 β 8.93701e-10 β
β 39.0 β 0.567143 β 5.06855e-10 β
β 40.0 β 0.567143 β 2.87461e-10 β
β 41.0 β 0.567143 β 1.6303e-10 β
β 42.0 β 0.567143 β 9.24626e-11 β
β 43.0 β 0.567143 β 5.24385e-11 β
ββββββββ΄ββββββββββββββ΄ββββββββββββββ
= - 1 .;println ("Iteration xβββ = g(xβ) with xβ = $ x1" )printErrors ( simple_iteration ( g, x1 ) )
Iteration xβββ = g(xβ) with xβ = -1.0
ββββββββ¬βββββββββββ¬ββββββββββββββ
β n β x[n] β |x[n] - Ξ©| β
ββββββββΌβββββββββββΌββββββββββββββ€
β 1.0 β -1.0 β 1.56714 β
β 2.0 β 2.71828 β 2.15114 β
β 3.0 β 0.065988 β 0.501155 β
β 4.0 β 0.936142 β 0.368999 β
β 5.0 β 0.392138 β 0.175006 β
β 6.0 β 0.675611 β 0.108468 β
β 7.0 β 0.508845 β 0.0582979 β
β 8.0 β 0.601189 β 0.034046 β
β 9.0 β 0.548159 β 0.018984 β
β 10.0 β 0.578013 β 0.0108695 β
β 11.0 β 0.561012 β 0.00613117 β
β 12.0 β 0.570631 β 0.00348793 β
β 13.0 β 0.565169 β 0.00197471 β
β 14.0 β 0.568264 β 0.00112105 β
β 15.0 β 0.566508 β 0.000635441 β
β 16.0 β 0.567504 β 0.0003605 β
β 17.0 β 0.566939 β 0.000204419 β
β 18.0 β 0.567259 β 0.000115946 β
β 19.0 β 0.567078 β 6.57544e-5 β
β 20.0 β 0.567181 β 3.72934e-5 β
β 21.0 β 0.567122 β 2.11503e-5 β
β 22.0 β 0.567155 β 1.19954e-5 β
β 23.0 β 0.567136 β 6.80306e-6 β
β 24.0 β 0.567147 β 3.85832e-6 β
β 25.0 β 0.567141 β 2.18822e-6 β
β 26.0 β 0.567145 β 1.24103e-6 β
β 27.0 β 0.567143 β 7.03844e-7 β
β 28.0 β 0.567144 β 3.99181e-7 β
β 29.0 β 0.567143 β 2.26393e-7 β
β 30.0 β 0.567143 β 1.28397e-7 β
β 31.0 β 0.567143 β 7.28195e-8 β
β 32.0 β 0.567143 β 4.12991e-8 β
β 33.0 β 0.567143 β 2.34225e-8 β
β 34.0 β 0.567143 β 1.32839e-8 β
β 35.0 β 0.567143 β 7.53388e-9 β
β 36.0 β 0.567143 β 4.27279e-9 β
β 37.0 β 0.567143 β 2.42329e-9 β
β 38.0 β 0.567143 β 1.37435e-9 β
β 39.0 β 0.567143 β 7.79454e-10 β
β 40.0 β 0.567143 β 4.42061e-10 β
β 41.0 β 0.567143 β 2.50713e-10 β
β 42.0 β 0.567143 β 1.42189e-10 β
β 43.0 β 0.567143 β 8.06427e-11 β
β 44.0 β 0.567143 β 4.5735e-11 β
ββββββββ΄βββββββββββ΄ββββββββββββββ
= 2 .;println ("Iteration xβββ = h(xβ) with xβ = $ x1" )printErrors ( simple_iteration ( h, x1 ) )
Iteration xβββ = h(xβ) with xβ = 2.0
βββββββββ¬βββββββββββββββ¬βββββββββββββ
β n β x[n] β |x[n] - Ξ©| β
βββββββββΌβββββββββββββββΌβββββββββββββ€
β 1.0 β 2.0 β 1.43286 β
β 2.0 β -0.693147 β 1.26029 β
β 3.0 β 0.366513 β 0.20063 β
β 4.0 β 1.00372 β 0.436578 β
β 5.0 β -0.0037146 β 0.570858 β
β 6.0 β 5.59549 β 5.02834 β
β 7.0 β -1.72196 β 2.2891 β
β 8.0 β -0.543463 β 1.11061 β
β 9.0 β 0.609793 β 0.04265 β
β 10.0 β 0.494635 β 0.072508 β
β 11.0 β 0.703935 β 0.136791 β
β 12.0 β 0.35107 β 0.216073 β
β 13.0 β 1.04677 β 0.479627 β
β 14.0 β -0.0457093 β 0.612853 β
β 15.0 β 3.08545 β 2.51831 β
β 16.0 β -1.1267 β 1.69384 β
β 17.0 β -0.119292 β 0.686435 β
β 18.0 β 2.12618 β 1.55904 β
β 19.0 β -0.754329 β 1.32147 β
β 20.0 β 0.281927 β 0.285216 β
β 21.0 β 1.26611 β 0.698963 β
β 22.0 β -0.235946 β 0.803089 β
β 23.0 β 1.44415 β 0.877009 β
β 24.0 β -0.367522 β 0.934666 β
β 25.0 β 1.00097 β 0.433828 β
β 26.0 β -0.000970968 β 0.568114 β
β 27.0 β 6.93722 β 6.37007 β
β 28.0 β -1.9369 β 2.50404 β
β 29.0 β -0.661089 β 1.22823 β
β 30.0 β 0.413867 β 0.153277 β
β 31.0 β 0.882212 β 0.315068 β
β 32.0 β 0.125323 β 0.44182 β
β 33.0 β 2.07686 β 1.50972 β
β 34.0 β -0.730856 β 1.298 β
β 35.0 β 0.313538 β 0.253605 β
β 36.0 β 1.15983 β 0.59269 β
β 37.0 β -0.148277 β 0.71542 β
β 38.0 β 1.90868 β 1.34153 β
β 39.0 β -0.646409 β 1.21355 β
β 40.0 β 0.436322 β 0.130821 β
β 41.0 β 0.829374 β 0.262231 β
β 42.0 β 0.187084 β 0.380059 β
β 43.0 β 1.6762 β 1.10905 β
β 44.0 β -0.516528 β 1.08367 β
β 45.0 β 0.660626 β 0.0934825 β
β 46.0 β 0.414568 β 0.152576 β
β 47.0 β 0.880519 β 0.313376 β
β 48.0 β 0.127244 β 0.4399 β
β 49.0 β 2.06165 β 1.49451 β
β 50.0 β -0.723507 β 1.29065 β
β 51.0 β 0.323645 β 0.243499 β
β 52.0 β 1.12811 β 0.560966 β
β 53.0 β -0.120543 β 0.687686 β
β 54.0 β 2.11575 β 1.54861 β
β 55.0 β -0.749409 β 1.31655 β
β 56.0 β 0.28847 β 0.278673 β
β 57.0 β 1.24316 β 0.67602 β
β 58.0 β -0.217659 β 0.784802 β
β 59.0 β 1.52483 β 0.957682 β
β 60.0 β -0.42188 β 0.989023 β
β 61.0 β 0.863034 β 0.295891 β
β 62.0 β 0.147301 β 0.419842 β
β 63.0 β 1.91528 β 1.34813 β
β 64.0 β -0.649862 β 1.21701 β
β 65.0 β 0.430995 β 0.136148 β
β 66.0 β 0.841659 β 0.274516 β
β 67.0 β 0.17238 β 0.394763 β
β 68.0 β 1.75805 β 1.19091 β
β 69.0 β -0.564207 β 1.13135 β
β 70.0 β 0.572333 β 0.00519021 β
β 71.0 β 0.558033 β 0.00910987 β
β 72.0 β 0.583336 β 0.0161931 β
β 73.0 β 0.538991 β 0.0281521 β
β 74.0 β 0.618056 β 0.0509127 β
β 75.0 β 0.481176 β 0.0859671 β
β 76.0 β 0.731522 β 0.164379 β
β 77.0 β 0.312628 β 0.254515 β
β 78.0 β 1.16274 β 0.595597 β
β 79.0 β -0.15078 β 0.717923 β
β 80.0 β 1.89193 β 1.32479 β
β 81.0 β -0.6376 β 1.20474 β
β 82.0 β 0.450045 β 0.117099 β
β 83.0 β 0.798408 β 0.231265 β
β 84.0 β 0.225135 β 0.342008 β
β 85.0 β 1.49105 β 0.923911 β
β 86.0 β -0.399483 β 0.966627 β
β 87.0 β 0.917583 β 0.35044 β
β 88.0 β 0.0860119 β 0.481131 β
β 89.0 β 2.45327 β 1.88613 β
β 90.0 β -0.897422 β 1.46457 β
β 91.0 β 0.108229 β 0.458914 β
β 92.0 β 2.2235 β 1.65636 β
β 93.0 β -0.799084 β 1.36623 β
β 94.0 β 0.224289 β 0.342854 β
β 95.0 β 1.49482 β 0.927676 β
β 96.0 β -0.402006 β 0.969149 β
β 97.0 β 0.911289 β 0.344146 β
β 98.0 β 0.0928949 β 0.474248 β
β 99.0 β 2.37629 β 1.80914 β
β 100.0 β -0.865539 β 1.43268 β
βββββββββ΄βββββββββββββββ΄βββββββββββββ
= Ξ© + 1e-10 ;println ("Iteration xβββ = h(xβ) with xβ = $ x1" )printErrors ( simple_iteration ( h, x1 ) )
Iteration xβββ = h(xβ) with xβ = 0.5671432905097845
βββββββββ¬βββββββββββββ¬ββββββββββββββ
β n β x[n] β |x[n] - Ξ©| β
βββββββββΌβββββββββββββΌββββββββββββββ€
β 1.0 β 0.567143 β 1.0e-10 β
β 2.0 β 0.567143 β 1.76324e-10 β
β 3.0 β 0.567143 β 3.10897e-10 β
β 4.0 β 0.567143 β 5.48182e-10 β
β 5.0 β 0.567143 β 9.66566e-10 β
β 6.0 β 0.567143 β 1.70427e-9 β
β 7.0 β 0.567143 β 3.00501e-9 β
β 8.0 β 0.567143 β 5.2985e-9 β
β 9.0 β 0.567143 β 9.34244e-9 β
β 10.0 β 0.567143 β 1.64728e-8 β
β 11.0 β 0.567143 β 2.90452e-8 β
β 12.0 β 0.567143 β 5.12132e-8 β
β 13.0 β 0.567143 β 9.03003e-8 β
β 14.0 β 0.567143 β 1.5922e-7 β
β 15.0 β 0.567144 β 2.8074e-7 β
β 16.0 β 0.567143 β 4.95006e-7 β
β 17.0 β 0.567144 β 8.72807e-7 β
β 18.0 β 0.567142 β 1.53895e-6 β
β 19.0 β 0.567146 β 2.71352e-6 β
β 20.0 β 0.567139 β 4.78453e-6 β
β 21.0 β 0.567152 β 8.43622e-6 β
β 22.0 β 0.567128 β 1.48748e-5 β
β 23.0 β 0.56717 β 2.6228e-5 β
β 24.0 β 0.567097 β 4.62447e-5 β
β 25.0 β 0.567225 β 8.15431e-5 β
β 26.0 β 0.567 β 0.000143768 β
β 27.0 β 0.567397 β 0.000253528 β
β 28.0 β 0.566696 β 0.000446926 β
β 29.0 β 0.567932 β 0.00078834 β
β 30.0 β 0.565754 β 0.00138905 β
β 31.0 β 0.569596 β 0.00245222 β
β 32.0 β 0.562829 β 0.00431448 β
β 33.0 β 0.57478 β 0.00763648 β
β 34.0 β 0.553768 β 0.013375 β
β 35.0 β 0.591009 β 0.0238656 β
β 36.0 β 0.525924 β 0.0412191 β
β 37.0 β 0.642598 β 0.0754548 β
β 38.0 β 0.442236 β 0.124908 β
β 39.0 β 0.815912 β 0.248769 β
β 40.0 β 0.203449 β 0.363695 β
β 41.0 β 1.59234 β 1.0252 β
β 42.0 β -0.465206 β 1.03235 β
β 43.0 β 0.765276 β 0.198132 β
β 44.0 β 0.267519 β 0.299624 β
β 45.0 β 1.31856 β 0.751421 β
β 46.0 β -0.276543 β 0.843687 β
β 47.0 β 1.28539 β 0.718245 β
β 48.0 β -0.251061 β 0.818204 β
β 49.0 β 1.38206 β 0.814918 β
β 50.0 β -0.323576 β 0.890719 β
β 51.0 β 1.12832 β 0.561178 β
β 52.0 β -0.120731 β 0.687875 β
β 53.0 β 2.11419 β 1.54705 β
β 54.0 β -0.748671 β 1.31581 β
β 55.0 β 0.289456 β 0.277688 β
β 56.0 β 1.23975 β 0.67261 β
β 57.0 β -0.214912 β 0.782056 β
β 58.0 β 1.53753 β 0.970382 β
β 59.0 β -0.430174 β 0.997318 β
β 60.0 β 0.843565 β 0.276422 β
β 61.0 β 0.170118 β 0.397025 β
β 62.0 β 1.77126 β 1.20412 β
β 63.0 β -0.571691 β 1.13883 β
β 64.0 β 0.559156 β 0.00798729 β
β 65.0 β 0.581327 β 0.0141835 β
β 66.0 β 0.542442 β 0.0247011 β
β 67.0 β 0.611674 β 0.0445304 β
β 68.0 β 0.491556 β 0.0755869 β
β 69.0 β 0.710179 β 0.143035 β
β 70.0 β 0.342239 β 0.224905 β
β 71.0 β 1.07225 β 0.505104 β
β 72.0 β -0.0697564 β 0.6369 β
β 73.0 β 2.66275 β 2.0956 β
β 74.0 β -0.979358 β 1.5465 β
β 75.0 β 0.0208581 β 0.546285 β
β 76.0 β 3.87001 β 3.30287 β
β 77.0 β -1.35326 β 1.9204 β
β 78.0 β -0.302515 β 0.869658 β
β 79.0 β 1.19563 β 0.628482 β
β 80.0 β -0.178669 β 0.745812 β
β 81.0 β 1.72222 β 1.15508 β
β 82.0 β -0.543614 β 1.11076 β
β 83.0 β 0.609516 β 0.0423726 β
β 84.0 β 0.49509 β 0.072053 β
β 85.0 β 0.703015 β 0.135872 β
β 86.0 β 0.352377 β 0.214767 β
β 87.0 β 1.04305 β 0.475911 β
β 88.0 β -0.0421532 β 0.609296 β
β 89.0 β 3.16644 β 2.5993 β
β 90.0 β -1.15261 β 1.71975 β
β 91.0 β -0.142028 β 0.709172 β
β 92.0 β 1.95173 β 1.38459 β
β 93.0 β -0.668715 β 1.23586 β
β 94.0 β 0.402397 β 0.164746 β
β 95.0 β 0.910317 β 0.343173 β
β 96.0 β 0.0939629 β 0.47318 β
β 97.0 β 2.36486 β 1.79771 β
β 98.0 β -0.860717 β 1.42786 β
β 99.0 β 0.149989 β 0.417154 β
β 100.0 β 1.89719 β 1.33005 β
βββββββββ΄βββββββββββββ΄ββββββββββββββ
We have therefore seen that x_{n+1} = g(x_n) converges to \Omega for all x_1 \in \{ -1, 1, 100 \} and y_{n+1} = g(y_n) diverges for y_1 \in \{ 1, \Omega + 10^{-10} \} .
We will show that \Omega is a stable fixed point of g but an unstable fixed point of h .
β Show that g:[\frac{1}{e}, 1]\to[\frac{1}{e},1] with |g'| < 1 ,
Answer.
The function g(x) = e^{-x} is decreasing with g(1) = \frac1e and g(\tfrac1e) = e^{-\tfrac1e} . Therefore, g( [\frac{1}{e}, 1] ) \subset [ \frac1e, e^{-\tfrac1e} ] \subset [\frac1e, 1] .
Moreover, we have |g'(x)| = |-e^{-x}| \leq e^{-\tfrac1e} < 1 for all x \in [\frac1e, 1] .
This can also be seen in the following plot:
println ( "[g(1), g(1/e)] = [1/e, " , g (1 / exp (1 )), "] β [1/e, 1]" )plot ( g, 1 / exp (1 ), 1 , label= L"y = g(x)" )hline! ( [1 / exp (1 ), 1 ] , linestyle= [: dash : dash], primary= false )scatter! ( [1 / exp (1 ), 1 ], [g (1 / exp (1 )), g (1 )], markersize= 5 , primary= false )
[g(1), g(1/e)] = [1/e, 0.6922006275553464] β [1/e, 1]
β Explain why this means that the iteration x_{n+1} = g(x_n) converges for all x_1 \in [\frac{1}{e}, 1] ,
Answer.
Since |g'(x)| \leq e^{-\frac1e} < 1 for all x \in [\frac{1}{e}, 1] (which we showed in part 3.), g is a contraction on [\frac{1}{e}, 1] . By the contraction mapping theorem, there exists unique fixed point \xi = g(\xi) \in [\frac{1}{e}, 1] and the iteration x_{n+1} = g(x_n) converges at least linearly for all x_1 \in [\frac{1}{e}, 1] .
β Now fix x_1 \in \mathbb R . Show that the sequence (x_n) is eventually contained in [0,1] and that g : [0,1] \to [\frac{1}{e}, 1] . Explain why this means that (x_n) converges for all x_1 \in \mathbb R ,
Answer.
For any x_1 \in \mathbb R , we have x_2 = g(x_1) = e^{-x_1} \geq 0 . Then, either x_2 \in [0,1] or x_2 > 1 . In the later case, x_3 = g(x_2) = e^{-x_2} \in [0, \frac1e] \subset [0,1] . That is, either x_2 \in [0,1] or x_3 \in [0,1] . Since g is decreasing, we have g([0,1])\in [\frac1e, 1 ] and so either x_3 \in [\frac1e, 1 ] or x_4 \in [\frac1e, 1 ] . By the previous part, the iteration x_{n+1} = g(x_n) started in [\frac1e,1] converges.
β What is the order of convergence of this method?
Answer.
We showed in lectures that
\begin{align}
\lim_{n\to\infty}\frac
{|x_{n+1} - \Omega |}
{|x_n - \Omega|}
= \left| g'(\Omega) \right|.
\end{align}
In this particular case |g'(\Omega)| = e^{-\Omega} = g(\Omega) = \Omega > 0 . As a result, the convergence is linear.
π» Confirm this numerically.
Answer.
By computing
\begin{align}
\frac
{|x_{n+1}-\Omega|}
{|x_n - \Omega|}
\end{align}
we see that this approaches a value in (0,1) as n\to\infty :
println ( "theoretical value of ΞΌ = " , Ξ©)= simple_iteration ( g, 1 . ) orderOfConvergence ( x, Ξ©, Ξ±= 1 )
theoretical value of ΞΌ = 0.5671432904097845
βββββββββββββ¬βββββββββββ¬ββββββββββββββ¬ββββββββββ¬ββββββββββββββ
β iteration β sequence β abs. error β ratio β ΞΌ (Ξ± = 1.0) β
βββββββββββββΌβββββββββββΌββββββββββββββΌββββββββββΌββββββββββββββ€
β 1.0 β 1.0 β 0.432857 β NaN β NaN β
β 2.0 β 0.367879 β 0.199264 β 1.92647 β 0.460346 β
β 3.0 β 0.692201 β 0.125057 β 1.28879 β 0.627597 β
β 4.0 β 0.500474 β 0.0666698 β 1.30256 β 0.533114 β
β 5.0 β 0.606244 β 0.0391002 β 1.19705 β 0.586476 β
β 6.0 β 0.545396 β 0.0217475 β 1.18097 β 0.556199 β
β 7.0 β 0.579612 β 0.012469 β 1.1453 β 0.573355 β
β 8.0 β 0.560115 β 0.00702783 β 1.13077 β 0.563622 β
β 9.0 β 0.571143 β 0.00399982 β 1.11368 β 0.569141 β
β 10.0 β 0.564879 β 0.00226394 β 1.10308 β 0.566011 β
β 11.0 β 0.568429 β 0.00128543 β 1.09293 β 0.567786 β
β 12.0 β 0.566415 β 0.000728557 β 1.0853 β 0.566779 β
β 13.0 β 0.567557 β 0.000413347 β 1.07845 β 0.56735 β
β 14.0 β 0.566909 β 0.000234378 β 1.07282 β 0.567026 β
β 15.0 β 0.567276 β 0.000132942 β 1.06784 β 0.56721 β
β 16.0 β 0.567068 β 7.5392e-5 β 1.06355 β 0.567106 β
β 17.0 β 0.567186 β 4.27597e-5 β 1.05974 β 0.567165 β
β 18.0 β 0.567119 β 2.42504e-5 β 1.05638 β 0.567131 β
β 19.0 β 0.567157 β 1.37536e-5 β 1.05337 β 0.56715 β
β 20.0 β 0.567135 β 7.8002e-6 β 1.05066 β 0.567139 β
β 21.0 β 0.567148 β 4.42385e-6 β 1.04822 β 0.567146 β
β 22.0 β 0.567141 β 2.50895e-6 β 1.046 β 0.567142 β
β 23.0 β 0.567145 β 1.42294e-6 β 1.04398 β 0.567144 β
β 24.0 β 0.567142 β 8.07008e-7 β 1.04213 β 0.567143 β
β 25.0 β 0.567144 β 4.5769e-7 β 1.04042 β 0.567144 β
β 26.0 β 0.567143 β 2.59576e-7 β 1.03885 β 0.567143 β
β 27.0 β 0.567143 β 1.47217e-7 β 1.0374 β 0.567143 β
β 28.0 β 0.567143 β 8.34929e-8 β 1.03605 β 0.567143 β
β 29.0 β 0.567143 β 4.73524e-8 β 1.0348 β 0.567143 β
β 30.0 β 0.567143 β 2.68556e-8 β 1.03363 β 0.567143 β
β 31.0 β 0.567143 β 1.5231e-8 β 1.03253 β 0.567143 β
β 32.0 β 0.567143 β 8.63815e-9 β 1.03151 β 0.567143 β
β 33.0 β 0.567143 β 4.89907e-9 β 1.03055 β 0.567143 β
β 34.0 β 0.567143 β 2.77847e-9 β 1.02964 β 0.567144 β
β 35.0 β 0.567143 β 1.57579e-9 β 1.02879 β 0.567143 β
β 36.0 β 0.567143 β 8.93701e-10 β 1.02798 β 0.567144 β
β 37.0 β 0.567143 β 5.06855e-10 β 1.02722 β 0.567142 β
β 38.0 β 0.567143 β 2.87461e-10 β 1.0265 β 0.567145 β
β 39.0 β 0.567143 β 1.6303e-10 β 1.02581 β 0.56714 β
β 40.0 β 0.567143 β 9.24626e-11 β 1.02516 β 0.56715 β
β 41.0 β 0.567143 β 5.24385e-11 β 1.02455 β 0.567132 β
βββββββββββββ΄βββββββββββ΄ββββββββββββββ΄ββββββββββ΄ββββββββββββββ
β Show that (y_n) doesnβt converge to \Omega for any choice of y_1 \not= \Omega
Answer.
We know that \Omega = h(\Omega) . The derivative is h'(x) = \frac1{1/x} \left(-\frac{1}{x^2}\right) = - \frac{1}{x} . Therefore, we have \left|h'(\Omega)\right| = \frac{1}{\Omega} \approx 1.76 (see code below) and in a small interval around \Omega , we have |h'| > 1 . By a result from lectures, this means that the sequence y_{n+1} = h(y_n) does not converge to \Omega .
(Recall the proof of this statement: if y_n \in I\setminus\{\xi\} for some small interval I around \xi for which |h'|>1 on I , then |y_{n+1} - \xi| = |h(y_n) - \xi| = |h'(\eta)| |x_{n}-\xi| > |x_{n}-\xi| and so the sequence is eventually leaves I )
B. Bisection Method
Suppose that f \colon [a,b] \to \mathbb R is continuous with f(a) f(b) \leq 0 . By the change of sign theorem, there exists \xi \in [a,b] for which f(\xi) = 0 . The bisection method can be described by the following pseudocode:
INPUT: f : function
a, b : interval on which f is defined
tol : tolerance
N : maximum number of iterations
OUTPUT: x : such that |f(x)| < tol
for n in 1:N
c = midpoint of a and b
if ( |f(c)| < tol )
OUTPUT: approximate root c
else if ( f changes sign on [a,c] )
b = c;
else
a = c;
end
OUTPUT: "root not found"
π» Implement this method in the following cell
function bisection ( f, a, b; tol= 1e-10 , N= 100 )= [];for n in 1 : N= (a+ b)/ 2 = push! (x, c)if ( abs ( f ( c ) ) < tol )println ("converged in " , n, " iterations" )return x elseif ( f (a)*f (c) < 0 )= celse = c end end @warn "did not converge after $ N iterations" return x end = x -> cos ( 10 * x^ 2 ); = bisection ( f, 0.0 , 1.0 );= sqrt (pi / 20 )= [1 : length (x) x (@.abs (x- ΞΎ))]pretty_table ( = ["n" , "x[n]" , "|x[n] - ΞΎ|" ],plot ( f , 0 , 1 , legend= false , lw = 2 )scatter! ( [x[end ]], [f (x[end ])], color= "red" )
converged in 35 iterations
ββββββ¬βββββββββββ¬ββββββββββββββ
β n β x[n] β |x[n] - ΞΎ| β
ββββββΌβββββββββββΌββββββββββββββ€
β 1 β 0.5 β 0.103667 β
β 2 β 0.25 β 0.146333 β
β 3 β 0.375 β 0.0213327 β
β 4 β 0.4375 β 0.0411673 β
β 5 β 0.40625 β 0.00991727 β
β 6 β 0.390625 β 0.00570773 β
β 7 β 0.398438 β 0.00210477 β
β 8 β 0.394531 β 0.00180148 β
β 9 β 0.396484 β 0.000151645 β
β 10 β 0.395508 β 0.000824917 β
β 11 β 0.395996 β 0.000336636 β
β 12 β 0.39624 β 9.24954e-5 β
β 13 β 0.396362 β 2.95749e-5 β
β 14 β 0.396301 β 3.14602e-5 β
β 15 β 0.396332 β 9.42651e-7 β
β 16 β 0.396347 β 1.43161e-5 β
β 17 β 0.396339 β 6.68674e-6 β
β 18 β 0.396336 β 2.87205e-6 β
β 19 β 0.396334 β 9.64697e-7 β
β 20 β 0.396333 β 1.10231e-8 β
β 21 β 0.396332 β 4.65814e-7 β
β 22 β 0.396333 β 2.27395e-7 β
β 23 β 0.396333 β 1.08186e-7 β
β 24 β 0.396333 β 4.85816e-8 β
β 25 β 0.396333 β 1.87792e-8 β
β 26 β 0.396333 β 3.87807e-9 β
β 27 β 0.396333 β 3.57251e-9 β
β 28 β 0.396333 β 1.52781e-10 β
β 29 β 0.396333 β 1.70986e-9 β
β 30 β 0.396333 β 7.78542e-10 β
β 31 β 0.396333 β 3.12881e-10 β
β 32 β 0.396333 β 8.00501e-11 β
β 33 β 0.396333 β 3.63652e-11 β
β 34 β 0.396333 β 2.18424e-11 β
β 35 β 0.396333 β 7.26141e-12 β
ββββββ΄βββββββββββ΄ββββββββββββββ
β Suppose there exists a unique root \xi \in [a,b] . What is the order of convergence of the bisection method? What is the asymptotic error constant?
Answer.
We have x_1 = \frac{a + b}{2} and thus |x_1 - \xi| \leq \frac{1}{2}(b-a) . Moreover, since the size of the interval containing \xi is halved at each iteration, we have |x_n - \xi| \leq \big(\frac{1}{2}\big)^{n} (b-a) . As a result, the bisection method converges linearly with asymptotic error constant \frac12 .
βπ» Plot f(x) = (1+x^2) \mathrm{sign}(x - \xi) and explain what happens when you apply the bisection method to f ,
Answer.
The algorithm should run as normal and converge to \xi : the same bracketing argument applies to f because it only has one point of discontinuity and a unique point at which the sign of f changes.
= 1 / Ο , x -> (1 + x^ 2 ) * sign (x - ΞΎ)= bisection ( f, 0.0 , 1.0 );= [1 : length (x) x (@.abs (x- ΞΎ))]pretty_table ( = ["n" , "x[n]" , "|x[n] - ΞΎ|" ],plot ( f , 0 , 1 , legend= false , lw = 2 )scatter! ( [x[end ]], [0 ], color= "red" )
converged in 54 iterations
ββββββ¬βββββββββββ¬ββββββββββββββ
β n β x[n] β |x[n] - ΞΎ| β
ββββββΌβββββββββββΌββββββββββββββ€
β 1 β 0.5 β 0.18169 β
β 2 β 0.25 β 0.0683099 β
β 3 β 0.375 β 0.0566901 β
β 4 β 0.3125 β 0.00580989 β
β 5 β 0.34375 β 0.0254401 β
β 6 β 0.328125 β 0.00981511 β
β 7 β 0.320312 β 0.00200261 β
β 8 β 0.316406 β 0.00190364 β
β 9 β 0.318359 β 4.94888e-5 β
β 10 β 0.317383 β 0.000927074 β
β 11 β 0.317871 β 0.000438792 β
β 12 β 0.318115 β 0.000194652 β
β 13 β 0.318237 β 7.25815e-5 β
β 14 β 0.318298 β 1.15463e-5 β
β 15 β 0.318329 β 1.89712e-5 β
β 16 β 0.318314 β 3.71245e-6 β
β 17 β 0.318306 β 3.91695e-6 β
β 18 β 0.31831 β 1.02248e-7 β
β 19 β 0.318312 β 1.8051e-6 β
β 20 β 0.318311 β 8.51426e-7 β
β 21 β 0.31831 β 3.74589e-7 β
β 22 β 0.31831 β 1.3617e-7 β
β 23 β 0.31831 β 1.6961e-8 β
β 24 β 0.31831 β 4.26436e-8 β
β 25 β 0.31831 β 1.28413e-8 β
β 26 β 0.31831 β 2.05988e-9 β
β 27 β 0.31831 β 5.3907e-9 β
β 28 β 0.31831 β 1.66541e-9 β
β 29 β 0.31831 β 1.97239e-10 β
β 30 β 0.31831 β 7.34083e-10 β
β 31 β 0.31831 β 2.68422e-10 β
β 32 β 0.31831 β 3.55913e-11 β
β 33 β 0.31831 β 8.08241e-11 β
β 34 β 0.31831 β 2.26164e-11 β
β 35 β 0.31831 β 6.48742e-12 β
β 36 β 0.31831 β 8.06449e-12 β
β 37 β 0.31831 β 7.88536e-13 β
β 38 β 0.31831 β 2.84944e-12 β
β 39 β 0.31831 β 1.03045e-12 β
β 40 β 0.31831 β 1.20959e-13 β
β 41 β 0.31831 β 3.33789e-13 β
β 42 β 0.31831 β 1.06415e-13 β
β 43 β 0.31831 β 7.27196e-15 β
β 44 β 0.31831 β 4.95715e-14 β
β 45 β 0.31831 β 2.11497e-14 β
β 46 β 0.31831 β 6.93889e-15 β
β 47 β 0.31831 β 1.66533e-16 β
β 48 β 0.31831 β 3.38618e-15 β
β 49 β 0.31831 β 1.60982e-15 β
β 50 β 0.31831 β 7.21645e-16 β
β 51 β 0.31831 β 2.77556e-16 β
β 52 β 0.31831 β 5.55112e-17 β
β 53 β 0.31831 β 5.55112e-17 β
β 54 β 0.31831 β 0.0 β
ββββββ΄βββββββββββ΄ββββββββββββββ
π» Use the bisection method to approximate \Omega = W(1) ,
= Ξ©, x -> x*exp (x)- 1 = bisection ( f, 0.0 , 1.0 );= [1 : length (x) x (@.abs (x- ΞΎ))]pretty_table ( = ["n" , "x[n]" , "|x[n] - ΞΎ|" ],plot ( f , 0 , 1 , legend= false , lw = 2 )scatter! ( [x[end ]], [0 ], color= "red" )
converged in 33 iterations
ββββββ¬βββββββββββ¬ββββββββββββββ
β n β x[n] β |x[n] - ΞΎ| β
ββββββΌβββββββββββΌββββββββββββββ€
β 1 β 0.5 β 0.0671433 β
β 2 β 0.75 β 0.182857 β
β 3 β 0.625 β 0.0578567 β
β 4 β 0.5625 β 0.00464329 β
β 5 β 0.59375 β 0.0266067 β
β 6 β 0.578125 β 0.0109817 β
β 7 β 0.570312 β 0.00316921 β
β 8 β 0.566406 β 0.00073704 β
β 9 β 0.568359 β 0.00121608 β
β 10 β 0.567383 β 0.000239522 β
β 11 β 0.566895 β 0.000248759 β
β 12 β 0.567139 β 4.61853e-6 β
β 13 β 0.567261 β 0.000117452 β
β 14 β 0.5672 β 5.64166e-5 β
β 15 β 0.567169 β 2.5899e-5 β
β 16 β 0.567154 β 1.06403e-5 β
β 17 β 0.567146 β 3.01086e-6 β
β 18 β 0.567142 β 8.03838e-7 β
β 19 β 0.567144 β 1.10351e-6 β
β 20 β 0.567143 β 1.49837e-7 β
β 21 β 0.567143 β 3.27e-7 β
β 22 β 0.567143 β 8.85818e-8 β
β 23 β 0.567143 β 3.06275e-8 β
β 24 β 0.567143 β 2.89771e-8 β
β 25 β 0.567143 β 8.25186e-10 β
β 26 β 0.567143 β 1.4076e-8 β
β 27 β 0.567143 β 6.6254e-9 β
β 28 β 0.567143 β 2.9001e-9 β
β 29 β 0.567143 β 1.03746e-9 β
β 30 β 0.567143 β 1.06137e-10 β
β 31 β 0.567143 β 3.59524e-10 β
β 32 β 0.567143 β 1.26694e-10 β
β 33 β 0.567143 β 1.02783e-11 β
ββββββ΄βββββββββββ΄ββββββββββββββ
π» Use the bisection method to approximate \sqrt{2} .
= sqrt (2 ), x -> x^ 2 - 2 = bisection ( f, 0.0 , 2.0 );= [1 : length (x) x (@.abs (x- ΞΎ))]pretty_table ( = ["n" , "x[n]" , "|x[n] - ΞΎ|" ],plot ( f , 0 , 2 , legend= false , lw = 2 )scatter! ( [x[end ]], [0 ], color= "red" )
converged in 30 iterations
ββββββ¬ββββββββββ¬ββββββββββββββ
β n β x[n] β |x[n] - ΞΎ| β
ββββββΌββββββββββΌββββββββββββββ€
β 1 β 1.0 β 0.414214 β
β 2 β 1.5 β 0.0857864 β
β 3 β 1.25 β 0.164214 β
β 4 β 1.375 β 0.0392136 β
β 5 β 1.4375 β 0.0232864 β
β 6 β 1.40625 β 0.00796356 β
β 7 β 1.42188 β 0.00766144 β
β 8 β 1.41406 β 0.000151062 β
β 9 β 1.41797 β 0.00375519 β
β 10 β 1.41602 β 0.00180206 β
β 11 β 1.41504 β 0.0008255 β
β 12 β 1.41455 β 0.000337219 β
β 13 β 1.41431 β 9.30783e-5 β
β 14 β 1.41418 β 2.89921e-5 β
β 15 β 1.41425 β 3.20431e-5 β
β 16 β 1.41422 β 1.52552e-6 β
β 17 β 1.4142 β 1.37333e-5 β
β 18 β 1.41421 β 6.10388e-6 β
β 19 β 1.41421 β 2.28918e-6 β
β 20 β 1.41421 β 3.81831e-7 β
β 21 β 1.41421 β 5.71843e-7 β
β 22 β 1.41421 β 9.50061e-8 β
β 23 β 1.41421 β 1.43413e-7 β
β 24 β 1.41421 β 2.42032e-8 β
β 25 β 1.41421 β 3.54014e-8 β
β 26 β 1.41421 β 5.59909e-9 β
β 27 β 1.41421 β 9.30207e-9 β
β 28 β 1.41421 β 1.85149e-9 β
β 29 β 1.41421 β 1.8738e-9 β
β 30 β 1.41421 β 1.11526e-11 β
ββββββ΄ββββββββββ΄ββββββββββββββ
