{"id":765,"date":"2013-09-16T12:00:42","date_gmt":"2013-09-16T17:00:42","guid":{"rendered":"https:\/\/blogs.mathworks.com\/cleve\/?p=765"},"modified":"2013-09-13T17:26:11","modified_gmt":"2013-09-13T22:26:11","slug":"tower-of-powers","status":"publish","type":"post","link":"https:\/\/blogs.mathworks.com\/cleve\/2013\/09\/16\/tower-of-powers\/","title":{"rendered":"Tower of Powers"},"content":{"rendered":"<div class=\"content\"><!--introduction--><p>I first investigated the tower of powers<\/p><p>$$ z^{z^{z^z}} $$<\/p><p>when I was in junior high school and was about 14 or 15 years old. I return to it every few decades.  This is another one of those episodes.<\/p><!--\/introduction--><h3>Contents<\/h3><div><ul><li><a href=\"#431a7cad-47d7-43c1-afb1-3503533dc2c6\">Associativity<\/a><\/li><li><a href=\"#77b1c4a3-8c5e-43fd-b03c-bec2f637e670\">Mathematical convention<\/a><\/li><li><a href=\"#5c01c22b-663d-45cc-a84a-7a21bc84291b\">WS Algorithm<\/a><\/li><li><a href=\"#11feff3b-c1c3-44c3-ae63-179add4819e1\">Two experiments<\/a><\/li><li><a href=\"#9adc78c5-e519-4aa4-8f1d-78443e9cbfd9\">sqrt(2)<\/a><\/li><li><a href=\"#6cd19342-2694-49f7-861f-165cf3d0c863\">Highest tower<\/a><\/li><li><a href=\"#f3b8a3c5-9ae3-4cb6-82ba-b72cd382e32b\">Analysis<\/a><\/li><li><a href=\"#bbb2051b-b557-46c9-8a44-56b303637ab7\">Stable Towers<\/a><\/li><li><a href=\"#70b3351c-be6a-453d-9a83-79d11295d931\">Towers Topple<\/a><\/li><li><a href=\"#e71ac1d6-ee51-4288-8cd2-5ed58827696a\">Small z<\/a><\/li><li><a href=\"#4c00875a-4cd5-4958-9aee-ffaa00c015f0\">Enter Lambert W<\/a><\/li><li><a href=\"#83ac6c1a-7a14-455c-aa09-c07cd8a2255c\">The Tower function<\/a><\/li><li><a href=\"#c9f3fb4a-0de1-4c38-ad44-949aae52456d\">Plot<\/a><\/li><li><a href=\"#fb8673a5-dd2a-4efe-9aad-b426712f455a\">Unstable<\/a><\/li><\/ul><\/div><h4>Associativity<a name=\"431a7cad-47d7-43c1-afb1-3503533dc2c6\"><\/a><\/h4><p>Does<\/p><p>$$ 2^{3^4} $$<\/p><p>mean<\/p><p>$$ {(2^3)}^4 = 8^4 = 2^{12},  \\ \\mbox{or}, \\ 2^{(3^4)} = 2^{81} $$<\/p><p>It makes a big difference.<\/p><p>Addition and multiplication are associative.  By that we mean<\/p><pre class=\"codeinput\">   2*3*4\r\n<\/pre><pre class=\"codeoutput\">ans =\r\n    24\r\n<\/pre><p>could have been computed from either<\/p><pre class=\"codeinput\">   (2*3)*4\r\n<\/pre><pre class=\"codeoutput\">ans =\r\n    24\r\n<\/pre><p>or<\/p><pre class=\"codeinput\">   2*(3*4)\r\n<\/pre><pre class=\"codeoutput\">ans =\r\n    24\r\n<\/pre><p>They're both the same.<\/p><p>But subtraction, division, and exponentiation are not associative. MATLAB associates from left to right.  So<\/p><pre class=\"codeinput\">   format <span class=\"string\">rat<\/span>\r\n   2\/3\/4\r\n<\/pre><pre class=\"codeoutput\">ans =\r\n       1\/6     \r\n<\/pre><p>was computed by<\/p><pre class=\"codeinput\">   (2\/3)\/4\r\n<\/pre><pre class=\"codeoutput\">ans =\r\n       1\/6     \r\n<\/pre><p>and not by<\/p><pre class=\"codeinput\">   2\/(3\/4)\r\n<\/pre><pre class=\"codeoutput\">ans =\r\n       8\/3     \r\n<\/pre><h4>Mathematical convention<a name=\"77b1c4a3-8c5e-43fd-b03c-bec2f637e670\"><\/a><\/h4><p>The usual mathematical convention for exponentiation is right associativity.<\/p><p>$$ a^{b^c} = a^{(b^c)} $$<\/p><p>This is because left associativity can be replaced by multiplication.<\/p><p>$$ {(a^b)}^c = a^{bc} $$<\/p><p>So, my tower of powers is<\/p><p>$$ z^{z^{z^z}} = z^{(z^{(z^z)})} $$<\/p><p>And, I want to make this tower infinitely high.<\/p><h4>WS Algorithm<a name=\"5c01c22b-663d-45cc-a84a-7a21bc84291b\"><\/a><\/h4><p>I've been stretching my typographical luck with superscripts and parentheses here. Let's settle down.<\/p><p>For a given value of $z$, let $y$ denote the infinite extension of the tower at the end of the last paragraph.  Then, since it's infinite, and associated to the right, we have<\/p><p>$$ y = z^y $$<\/p><p>Now that's at least easier to typeset.<\/p><p>There is an obvious way to try to compute it. I call this the <i>WS<\/i> or <i>World's Simplest<\/i> algorithm.<\/p><p>$$ y_1 = 1 $$<\/p><p>$$ \\mbox{for} \\ n = 1, ..., \\ \\ y_{n+1} = z^{y_n} $$<\/p><p>Of course, this is better known as the <i>Fixed Point<\/i> algorithm<\/p><p>$$ y_{n+1} = f(y_n) $$<\/p><p>for the function<\/p><p>$$ f(z) = z^y $$.<\/p><h4>Two experiments<a name=\"11feff3b-c1c3-44c3-ae63-179add4819e1\"><\/a><\/h4><p>Let's try two different values of $z$, $1.25$ and $1.50$.<\/p><pre class=\"codeinput\">   format <span class=\"string\">long<\/span>\r\n   z = 1.25\r\n   y = 1\r\n   <span class=\"keyword\">for<\/span> n = 1:30\r\n      y = z^y;\r\n      disp(y)\r\n   <span class=\"keyword\">end<\/span>\r\n   snapnow\r\n<\/pre><pre class=\"codeoutput\">z =\r\n   1.250000000000000\r\ny =\r\n     1\r\n   1.250000000000000\r\n   1.321714079300705\r\n   1.343034993578008\r\n   1.349439873733169\r\n   1.351369882459350\r\n   1.351952000918108\r\n   1.352127625454628\r\n   1.352180615675240\r\n   1.352196604529415\r\n   1.352201428918186\r\n   1.352202884606055\r\n   1.352203323838621\r\n   1.352203456370664\r\n   1.352203496360284\r\n   1.352203508426572\r\n   1.352203512067399\r\n   1.352203513165966\r\n   1.352203513497443\r\n   1.352203513597461\r\n   1.352203513627640\r\n   1.352203513636746\r\n   1.352203513639493\r\n   1.352203513640323\r\n   1.352203513640573\r\n   1.352203513640648\r\n   1.352203513640671\r\n   1.352203513640678\r\n   1.352203513640680\r\n   1.352203513640681\r\n   1.352203513640681\r\n<\/pre><pre class=\"codeinput\">   format <span class=\"string\">long<\/span>\r\n   z = 1.5\r\n   y = 1\r\n   <span class=\"keyword\">for<\/span> n = 1:15\r\n      y = z^y;\r\n      disp(y)\r\n   <span class=\"keyword\">end<\/span>\r\n   snapnow\r\n<\/pre><pre class=\"codeoutput\">z =\r\n   1.500000000000000\r\ny =\r\n     1\r\n   1.500000000000000\r\n   1.837117307087384\r\n   2.106203352148880\r\n   2.349005318611934\r\n   2.592025704907509\r\n   2.860441497460648\r\n   3.189324761899238\r\n   3.644283987904574\r\n   4.382546732246116\r\n   5.911914873330858\r\n  10.990982932689437\r\n  86.181891743310985\r\n     1.499263005586202e+15\r\n   Inf\r\n   Inf\r\n<\/pre><p>It converges for  $z = 1.25$, but not for $z = 1.5$. Something important is happening between $1.25$ and $1.5$.<\/p><h4>sqrt(2)<a name=\"9adc78c5-e519-4aa4-8f1d-78443e9cbfd9\"><\/a><\/h4><p>Let's try $z = \\sqrt{2}$<\/p><pre class=\"codeinput\">   format <span class=\"string\">long<\/span>\r\n   z = sqrt(2)\r\n   y = 1\r\n   <span class=\"keyword\">for<\/span> n = 1:30\r\n      y = z^y;\r\n      disp(y)\r\n   <span class=\"keyword\">end<\/span>\r\n<\/pre><pre class=\"codeoutput\">z =\r\n   1.414213562373095\r\ny =\r\n     1\r\n   1.414213562373095\r\n   1.632526919438153\r\n   1.760839555880029\r\n   1.840910869291011\r\n   1.892712696828511\r\n   1.926999701847101\r\n   1.950034773805818\r\n   1.965664886517319\r\n   1.976341754409703\r\n   1.983668399303822\r\n   1.988711773413954\r\n   1.992190882947058\r\n   1.994594450712102\r\n   1.996256666265859\r\n   1.997407001141337\r\n   1.998203477508703\r\n   1.998755133084593\r\n   1.999137310119392\r\n   1.999402118324998\r\n   1.999585622935682\r\n   1.999712796329641\r\n   1.999800935492971\r\n   1.999862023757784\r\n   1.999904364443337\r\n   1.999933711582101\r\n   1.999954052897822\r\n   1.999968152149245\r\n   1.999977924873871\r\n   1.999984698747096\r\n   1.999989394007813\r\n<\/pre><p>It's converging to $y = 2$, because that's a fixed point.<\/p><p>$$ (\\sqrt{2})^2 = 2 $$<\/p><p>So the crucial point is between $1.4$ and $1.5$.<\/p><h4>Highest tower<a name=\"6cd19342-2694-49f7-861f-165cf3d0c863\"><\/a><\/h4><p>I remember being fascinating when I first learned that the largest possible value of $z$ is a value just a little larger than $\\sqrt{2}$. It is the surprising quantity<\/p><p>$$ \\bar{\\mbox{e}} = \\mbox{e}^{1\/\\mbox{e}} $$<\/p><pre class=\"codeinput\">   ebar = exp(1)^exp(-1)\r\n<\/pre><pre class=\"codeoutput\">ebar =\r\n   1.444667861009766\r\n<\/pre><p>This works because $\\bar{\\mbox{e}}$ is another fixed point.<\/p><p>$$\\bar{\\mbox{e}}^{\\mbox{e}}=(\\mbox{e}^{1\/{\\mbox{e}}})^{\\mbox{e}}=\\mbox{e}$$<\/p><p>Let's see why $\\bar{\\mbox{e}}$ is the largest possible value of $z$.<\/p><h4>Analysis<a name=\"f3b8a3c5-9ae3-4cb6-82ba-b72cd382e32b\"><\/a><\/h4><p>Analysis of the WS algorithm is provided by the mean value theorem.<\/p><p>$$ f(y_n) = f(y_{n-1}) + f'(\\xi_n)(y_n-y_{n-1}) $$<\/p><p>for some $\\xi_n$ between $y_n$ and $y_{n-1}$.<\/p><p>In our case, $f'(z) = \\ln{z}\\ z^y$, so<\/p><p>$$ y_{n+1} = y_n + (\\ln{z})(z^{\\xi_n})(y_n-y_{n-1})$$<\/p><p>This tells us that<\/p><p>$$ |y_{n+1}-y_n| \\le C |y_n-y_{n-1}|$$<\/p><p>for the quantity<\/p><p>$$C = \\mbox{max} |\\ln{z}||z^y|$$<\/p><p>where the maximum is taken over the $y$ s and $z$ s in a region including the iterates. If $C$ is less than 1, the iterates get closer together, and so the iteration is a <i>contraction<\/i> that converges.<\/p><p>Let's do a numerical search by making a contour plot of $(\\ln{z})(z^y)$. Use green for the region where $C$ is less than $1$, and red and blue for the regions where it is larger than 1.<\/p><pre class=\"codeinput\">   [Z,Y] = ndgrid(.01:.01:2.0,.01:.01:4.0);\r\n   R = log(Z).*Z.^Y;\r\n   contourf(Z,Y,R,[-5:4:-1 1:20:21]);\r\n   colormap([0 0 .75;0 .75 .25;1 0 0]);\r\n   e = exp(1);\r\n   line([ebar ebar],[e e],<span class=\"string\">'marker'<\/span>,<span class=\"string\">'.'<\/span>,<span class=\"string\">'markersize'<\/span>,24,<span class=\"string\">'color'<\/span>,<span class=\"string\">'yellow'<\/span>)\r\n   set(gca,<span class=\"string\">'xtick'<\/span>,[.01 1\/e 1 e^(1\/e) 2])\r\n   set(gca,<span class=\"string\">'xticklabel'<\/span>,{<span class=\"string\">'0'<\/span>,<span class=\"string\">'1\/e'<\/span>,<span class=\"string\">'1'<\/span>,<span class=\"string\">'e^(1\/e)'<\/span>,<span class=\"string\">'2'<\/span>})\r\n   set(gca,<span class=\"string\">'ytick'<\/span>,[.01 1 2 e 4])\r\n   set(gca,<span class=\"string\">'yticklabel'<\/span>,{<span class=\"string\">'0'<\/span>,<span class=\"string\">'1'<\/span>,<span class=\"string\">'2'<\/span>,<span class=\"string\">'e'<\/span>,<span class=\"string\">'4'<\/span>})\r\n   xlabel(<span class=\"string\">'z'<\/span>), ylabel(<span class=\"string\">'y'<\/span>)\r\n   axis <span class=\"string\">square<\/span>\r\n<\/pre><img decoding=\"async\" vspace=\"5\" hspace=\"5\" src=\"https:\/\/blogs.mathworks.com\/images\/cleve\/tower_01.png\" alt=\"\"> <h4>Stable Towers<a name=\"bbb2051b-b557-46c9-8a44-56b303637ab7\"><\/a><\/h4><p>If we stay in the green region, where $C$ is less than 1, the WS iteration converges, and the tower is stable.  That includes all values of $z$ in the interval<\/p><p>$$ 1\/\\mbox{e} \\le z \\le \\mbox{e}^{1\/\\mbox{e}} $$<\/p><h4>Towers Topple<a name=\"70b3351c-be6a-453d-9a83-79d11295d931\"><\/a><\/h4><p>The crucial value that I encountered in junior high is the yellow dot, $z=\\bar{\\mbox{e}}=\\mbox{e}^{1\/\\mbox{e}}$ with $y=\\mbox{e}$. As we approach the dot from the stable region, $f'(z)$ approaches $f'(\\bar{\\mbox{e}})=1$ and the WS algorithm converges incredibly slowly.<\/p><p>To the right of $z=\\bar{\\mbox{e}}$, the iteration encounters values of $y$ in the red danger region where C is larger than 1, so it diverges rapidly and the tower collapses quickly, as we saw at $z=1.5$.<\/p><h4>Small z<a name=\"e71ac1d6-ee51-4288-8cd2-5ed58827696a\"><\/a><\/h4><p>Values of $z$ less than $1\/\\mbox{e}$, which is about $0.37$, are a little more delicate.  We might be in the blue region, or might not. It turns out that values of $z$ greater than $\\mbox{e}^{\\mbox{-e}}$, which is about $0.066$, still encounter values of $y$ larger enough to keep the iteration in the safe green region, so it converges.<\/p><p>But now try $z = 0.05$ for a few iterations.<\/p><pre class=\"codeinput\">   z = 0.05\r\n   y = 1\r\n   <span class=\"keyword\">for<\/span> k = 1:20, y = z^y; disp(y), <span class=\"keyword\">end<\/span>\r\n<\/pre><pre class=\"codeoutput\">z =\r\n   0.050000000000000\r\ny =\r\n     1\r\n   0.050000000000000\r\n   0.860891659331735\r\n   0.075849745545982\r\n   0.796741072475616\r\n   0.091921259400909\r\n   0.759290007184149\r\n   0.102834993569945\r\n   0.734866735033395\r\n   0.110641061786732\r\n   0.717881331845427\r\n   0.116416584607138\r\n   0.705567440560545\r\n   0.120791283462031\r\n   0.696381005849372\r\n   0.124161635028397\r\n   0.689385252311344\r\n   0.126791198828070\r\n   0.683975974940914\r\n   0.128862555680671\r\n   0.679744887323325\r\n<\/pre><p>It doesn't look like it wants to converge. Run 300 more iterations without any output.<\/p><pre class=\"codeinput\">   <span class=\"keyword\">for<\/span> k = 1:300, y = z^y; <span class=\"keyword\">end<\/span>\r\n<\/pre><p>Now check on how we're doing after a total of 320 iterations.<\/p><pre class=\"codeinput\">   <span class=\"keyword\">for<\/span> k = 1:10, y = z^y; disp(y), <span class=\"keyword\">end<\/span>\r\n<\/pre><pre class=\"codeoutput\">   0.137359395737956\r\n   0.662660838900549\r\n   0.137359395737956\r\n   0.662660838900549\r\n   0.137359395737956\r\n   0.662660838900549\r\n   0.137359395737956\r\n   0.662660838900549\r\n   0.137359395737956\r\n   0.662660838900549\r\n<\/pre><p>It has settled into a pattern, oscillating back and forth between two values.  This is a simple example of a <i>period doubling bifurcation<\/i>.<\/p><p>This is the behavior for all $z$ less than $\\mbox{e}^{\\mbox{-e}}$. The values of $y$ encountered put us in the blue region.  The fixed point iteration is not a contraction.  It settles into a stable two-cycle. The distance between the two limit points increases as $z$ gets smaller.<\/p><p>A tower based on a value of $z$ less than $\\mbox{e}^{\\mbox{-e}}$ will have its odd numbered floors one height and its even numbered floors another.<\/p><h4>Enter Lambert W<a name=\"4c00875a-4cd5-4958-9aee-ffaa00c015f0\"><\/a><\/h4><p>Here is another way to look at the tower of powers. It involves the Lambert W function from <a href=\"https:\/\/blogs.mathworks.com\/cleve\/the-lambert-w-function\">my previous blog<\/a>.<\/p><p>$$ y = z^y = \\mbox{e}^{y\\ln{z}} $$<\/p><p>$$ y \\mbox{e}^{-y\\ln{z}} = 1 $$<\/p><p>$$ -y \\ln{z} \\ \\mbox{e}^{-y\\ln{z}} = -\\ln{z} $$<\/p><p>$$ -y \\ln{z} = W(-\\ln{z}) $$<\/p><p>$$ y = \\frac{W(-\\ln{z})}{-\\ln{z}} $$<\/p><p>So, the Lambert W function provides a value of the infinite tower for any value of $z$.  I didn't know this when I was in junior high school.<\/p><h4>The Tower function<a name=\"83ac6c1a-7a14-455c-aa09-c07cd8a2255c\"><\/a><\/h4><p>Let's try this formulation for a few values of $z$, including some smaller than $\\mbox{e}^{\\mbox{-e}}$.<\/p><pre class=\"codeinput\">  T = @(z) Lambert_W(-log(z)).\/(-log(z))\r\n  z = [1.25 sqrt(2) ebar 1.5 e^(-e) 0.05 10.^(-(3:3:15))]';\r\n  disp(<span class=\"string\">'         z                    T(z)'<\/span>)\r\n  disp([z T(z)])\r\n<\/pre><pre class=\"codeoutput\">T = \r\n    @(z)Lambert_W(-log(z)).\/(-log(z))\r\n         z                    T(z)\r\n   1.250000000000000   1.352203513640681\r\n   1.414213562373095   2.000000000000000\r\n   1.444667861009766   2.718281787975775\r\n   1.500000000000000                 NaN\r\n   0.065988035845313   0.367879441171442\r\n   0.050000000000000   0.350224852743194\r\n   0.001000000000000   0.219513151627722\r\n   0.000001000000000   0.141526856553182\r\n   0.000000001000000   0.107583717152717\r\n   0.000000000001000   0.087971500166559\r\n   0.000000000000001   0.074997053328098\r\n<\/pre><h4>Plot<a name=\"c9f3fb4a-0de1-4c38-ad44-949aae52456d\"><\/a><\/h4><p>And now a plot.<\/p><pre class=\"codeinput\">  ezplot(T,[0 ebar])\r\n  axis([0 ebar 0 e])\r\n  title(<span class=\"string\">'Tower of Powers'<\/span>)\r\n<\/pre><img decoding=\"async\" vspace=\"5\" hspace=\"5\" src=\"https:\/\/blogs.mathworks.com\/images\/cleve\/tower_02.png\" alt=\"\"> <h4>Unstable<a name=\"fb8673a5-dd2a-4efe-9aad-b426712f455a\"><\/a><\/h4><p>The value of $y$ obtained from the Lambert W formulation is a fixed point of $f(z)$, but a possibly unstable one.  If we have a value of $z$ less than $\\mbox{e}^{\\mbox{-e}}$ and we start the fixed point algorithm at the value of $y$ given by the Lambert W formula exactly, it will stay there.  But if we start with any other value of $y$, the iterates will eventually find their way to the pair of bifurcation attractors.<\/p><script language=\"JavaScript\"> <!-- \r\n    function grabCode_e9a2d52dcf9a4d1c856da167477af7f1() {\r\n        \/\/ Remember the title so we can use it in the new page\r\n        title = document.title;\r\n\r\n        \/\/ Break up these strings so that their presence\r\n        \/\/ in the Javascript doesn't mess up the search for\r\n        \/\/ the MATLAB code.\r\n        t1='e9a2d52dcf9a4d1c856da167477af7f1 ' + '##### ' + 'SOURCE BEGIN' + ' #####';\r\n        t2='##### ' + 'SOURCE END' + ' #####' + ' e9a2d52dcf9a4d1c856da167477af7f1';\r\n    \r\n        b=document.getElementsByTagName('body')[0];\r\n        i1=b.innerHTML.indexOf(t1)+t1.length;\r\n        i2=b.innerHTML.indexOf(t2);\r\n \r\n        code_string = b.innerHTML.substring(i1, i2);\r\n        code_string = code_string.replace(\/REPLACE_WITH_DASH_DASH\/g,'--');\r\n\r\n        \/\/ Use \/x3C\/g instead of the less-than character to avoid errors \r\n        \/\/ in the XML parser.\r\n        \/\/ Use '\\x26#60;' instead of '<' so that the XML parser\r\n        \/\/ doesn't go ahead and substitute the less-than character. \r\n        code_string = code_string.replace(\/\\x3C\/g, '\\x26#60;');\r\n\r\n        copyright = 'Copyright 2013 The MathWorks, Inc.';\r\n\r\n        w = window.open();\r\n        d = w.document;\r\n        d.write('<pre>\\n');\r\n        d.write(code_string);\r\n\r\n        \/\/ Add copyright line at the bottom if specified.\r\n        if (copyright.length > 0) {\r\n            d.writeln('');\r\n            d.writeln('%%');\r\n            if (copyright.length > 0) {\r\n                d.writeln('% _' + copyright + '_');\r\n            }\r\n        }\r\n\r\n        d.write('<\/pre>\\n');\r\n\r\n        d.title = title + ' (MATLAB code)';\r\n        d.close();\r\n    }   \r\n     --> <\/script><p style=\"text-align: right; font-size: xx-small; font-weight:lighter;   font-style: italic; color: gray\"><br><a href=\"javascript:grabCode_e9a2d52dcf9a4d1c856da167477af7f1()\"><span style=\"font-size: x-small;        font-style: italic;\">Get \r\n      the MATLAB code <noscript>(requires JavaScript)<\/noscript><\/span><\/a><br><br>\r\n      Published with MATLAB&reg; R2013b<br><\/p><p class=\"footer\"><br>\r\n      Published with MATLAB&reg; R2013b<br><\/p><\/div><!--\r\ne9a2d52dcf9a4d1c856da167477af7f1 ##### SOURCE BEGIN #####\r\n%% Tower of Powers\r\n% I first investigated the tower of powers\r\n%\r\n% $$ z^{z^{z^z}} $$\r\n%\r\n% when I was in junior high school and was about 14 or 15 years old.\r\n% I return to it every few decades.  This is another one of those episodes.\r\n\r\n%% Associativity\r\n% Does\r\n%\r\n% $$ 2^{3^4} $$\r\n%\r\n% mean\r\n%\r\n% $$ {(2^3)}^4 = 8^4 = 2^{12},  \\ \\mbox{or}, \\ 2^{(3^4)} = 2^{81} $$\r\n%\r\n% It makes a big difference.\r\n%\r\n\r\n%%\r\n% Addition and multiplication are associative.  By that we mean\r\n%\r\n   2*3*4\r\n%%\r\n% could have been computed from either\r\n%\r\n   (2*3)*4\r\n%%\r\n% or\r\n%\r\n   2*(3*4)\r\n%%\r\n% They're both the same.\r\n%\r\n\r\n%%\r\n% But subtraction, division, and exponentiation are not associative.\r\n% MATLAB associates from left to right.  So\r\n%\r\n   format rat\r\n   2\/3\/4\r\n%%\r\n% was computed by\r\n%\r\n   (2\/3)\/4\r\n%%\r\n% and not by\r\n%\r\n   2\/(3\/4)\r\n\r\n%% Mathematical convention\r\n% The usual mathematical convention for exponentiation is right associativity.\r\n%\r\n% $$ a^{b^c} = a^{(b^c)} $$\r\n%\r\n% This is because left associativity can be replaced by multiplication.\r\n%\r\n% $$ {(a^b)}^c = a^{bc} $$\r\n\r\n%%\r\n% So, my tower of powers is\r\n%\r\n% $$ z^{z^{z^z}} = z^{(z^{(z^z)})} $$\r\n%\r\n\r\n%%\r\n% And, I want to make this tower infinitely high.\r\n\r\n%% WS Algorithm\r\n% I've been stretching my typographical luck with superscripts and parentheses here.\r\n% Let's settle down.\r\n\r\n%%\r\n% For a given value of $z$, let $y$ denote the infinite extension of the\r\n% tower at the end of the last paragraph.  Then, since it's infinite, and\r\n% associated to the right, we have\r\n%\r\n% $$ y = z^y $$\r\n%\r\n% Now that's at least easier to typeset.\r\n\r\n%%\r\n% There is an obvious way to try to compute it.\r\n% I call this the _WS_ or _World's Simplest_ algorithm.\r\n%\r\n% $$ y_1 = 1 $$\r\n%\r\n% $$ \\mbox{for} \\ n = 1, ..., \\ \\ y_{n+1} = z^{y_n} $$\r\n\r\n%%\r\n% Of course, this is better known as the _Fixed Point_ algorithm\r\n%\r\n% $$ y_{n+1} = f(y_n) $$\r\n%\r\n% for the function\r\n%\r\n% $$ f(z) = z^y $$.\r\n\r\n%% Two experiments\r\n% Let's try two different values of $z$, $1.25$ and $1.50$.\r\n\r\n   format long\r\n   z = 1.25\r\n   y = 1\r\n   for n = 1:30\r\n      y = z^y;\r\n      disp(y)\r\n   end\r\n   snapnow\r\n\r\n%%\r\n\r\n   format long\r\n   z = 1.5\r\n   y = 1\r\n   for n = 1:15\r\n      y = z^y;\r\n      disp(y)\r\n   end\r\n   snapnow\r\n\r\n%%\r\n% It converges for  $z = 1.25$, but not for $z = 1.5$.\r\n% Something important is happening between $1.25$ and $1.5$.\r\n\r\n%% sqrt(2)\r\n% Let's try $z = \\sqrt{2}$\r\n\r\n   format long\r\n   z = sqrt(2)\r\n   y = 1\r\n   for n = 1:30\r\n      y = z^y;\r\n      disp(y)\r\n   end\r\n\r\n%%\r\n% It's converging to $y = 2$, because that's a fixed point.\r\n%\r\n% $$ (\\sqrt{2})^2 = 2 $$\r\n%\r\n% So the crucial point is between $1.4$ and $1.5$.\r\n\r\n%% Highest tower\r\n% I remember being fascinating when I first learned that the largest\r\n% possible value of $z$ is a value just a little larger than $\\sqrt{2}$.\r\n% It is the surprising quantity\r\n%\r\n% $$ \\bar{\\mbox{e}} = \\mbox{e}^{1\/\\mbox{e}} $$\r\n%\r\n   ebar = exp(1)^exp(-1)\r\n\r\n%%\r\n% This works because $\\bar{\\mbox{e}}$ is another fixed point.\r\n%\r\n% $$\\bar{\\mbox{e}}^{\\mbox{e}}=(\\mbox{e}^{1\/{\\mbox{e}}})^{\\mbox{e}}=\\mbox{e}$$\r\n\r\n%%\r\n% Let's see why $\\bar{\\mbox{e}}$ is the largest possible value\r\n% of $z$.\r\n\r\n%% Analysis\r\n% Analysis of the WS algorithm is provided by the mean value theorem.\r\n%\r\n% $$ f(y_n) = f(y_{n-1}) + f'(\\xi_n)(y_n-y_{n-1}) $$\r\n%\r\n% for some $\\xi_n$ between $y_n$ and $y_{n-1}$.\r\n\r\n%%\r\n% In our case, $f'(z) = \\ln{z}\\ z^y$, so\r\n%\r\n% $$ y_{n+1} = y_n + (\\ln{z})(z^{\\xi_n})(y_n-y_{n-1})$$\r\n%\r\n% This tells us that\r\n%\r\n% $$ |y_{n+1}-y_n| \\le C |y_n-y_{n-1}|$$\r\n%\r\n% for the quantity\r\n% \r\n% $$C = \\mbox{max} |\\ln{z}||z^y|$$\r\n%\r\n% where the maximum is taken over the $y$ s and $z$ s in a region\r\n% including the iterates.\r\n% If $C$ is less than 1, the iterates get closer together, and so the iteration\r\n% is a _contraction_ that converges.\r\n\r\n%%\r\n% Let's do a numerical search by making a contour plot of $(\\ln{z})(z^y)$.\r\n% Use green for the region where $C$ is less than $1$, and red and blue\r\n% for the regions where it is larger than 1.\r\n\r\n   [Z,Y] = ndgrid(.01:.01:2.0,.01:.01:4.0); \r\n   R = log(Z).*Z.^Y;\r\n   contourf(Z,Y,R,[-5:4:-1 1:20:21]);\r\n   colormap([0 0 .75;0 .75 .25;1 0 0]);\r\n   e = exp(1);\r\n   line([ebar ebar],[e e],'marker','.','markersize',24,'color','yellow')\r\n   set(gca,'xtick',[.01 1\/e 1 e^(1\/e) 2])\r\n   set(gca,'xticklabel',{'0','1\/e','1','e^(1\/e)','2'})\r\n   set(gca,'ytick',[.01 1 2 e 4])\r\n   set(gca,'yticklabel',{'0','1','2','e','4'})\r\n   xlabel('z'), ylabel('y')\r\n   axis square\r\n\r\n%% Stable Towers\r\n% If we stay in the green region, where $C$ is less than 1, the WS iteration\r\n% converges, and the tower is stable.  That includes all values of $z$\r\n% in the interval\r\n%\r\n% $$ 1\/\\mbox{e} \\le z \\le \\mbox{e}^{1\/\\mbox{e}} $$\r\n\r\n%% Towers Topple\r\n% The crucial value that I encountered in junior high is the yellow dot,\r\n% $z=\\bar{\\mbox{e}}=\\mbox{e}^{1\/\\mbox{e}}$ with $y=\\mbox{e}$.\r\n% As we approach the dot from the stable region, $f'(z)$ approaches\r\n% $f'(\\bar{\\mbox{e}})=1$ and the WS algorithm converges\r\n% incredibly slowly.\r\n\r\n%%\r\n% To the right of $z=\\bar{\\mbox{e}}$, the iteration encounters values\r\n% of $y$ in the red danger region where C is larger than 1, so it diverges\r\n% rapidly and the tower collapses quickly, as we saw at $z=1.5$.\r\n\r\n%% Small z\r\n% Values of $z$ less than $1\/\\mbox{e}$, which is about $0.37$, are a little\r\n% more delicate.  We might be in the blue region, or might not.\r\n% It turns out that values of $z$ greater than $\\mbox{e}^{\\mbox{-e}}$,\r\n% which is about $0.066$, still encounter values of $y$ larger enough to keep\r\n% the iteration in the safe green region, so it converges.\r\n\r\n%%\r\n% But now try $z = 0.05$ for a few iterations.\r\n%\r\n   z = 0.05\r\n   y = 1\r\n   for k = 1:20, y = z^y; disp(y), end\r\n\r\n%%\r\n% It doesn't look like it wants to converge.\r\n% Run 300 more iterations without any output.\r\n\r\n   for k = 1:300, y = z^y; end\r\n\r\n%%\r\n% Now check on how we're doing after a total of 320 iterations.\r\n\r\n   for k = 1:10, y = z^y; disp(y), end\r\n\r\n%%\r\n% It has settled into a pattern, oscillating back and forth between\r\n% two values.  This is a simple example of a _period doubling bifurcation_.\r\n\r\n%% \r\n% This is the behavior for all $z$ less than $\\mbox{e}^{\\mbox{-e}}$.\r\n% The values of $y$ encountered put us in the blue region.  The fixed point\r\n% iteration is not a contraction.  It settles into a stable two-cycle.\r\n% The distance between the two limit points increases as $z$ gets smaller.\r\n\r\n%%\r\n% A tower based on a value of $z$ less than $\\mbox{e}^{\\mbox{-e}}$\r\n% will have its odd numbered floors one height and its even numbered floors\r\n% another.\r\n\r\n%% Enter Lambert W\r\n% Here is another way to look at the tower of powers.\r\n% It involves the Lambert W function from\r\n% <https:\/\/blogs.mathworks.com\/cleve\/the-lambert-w-function my previous blog>.\r\n%\r\n% $$ y = z^y = \\mbox{e}^{y\\ln{z}} $$\r\n%\r\n% $$ y \\mbox{e}^{-y\\ln{z}} = 1 $$\r\n%\r\n% $$ -y \\ln{z} \\ \\mbox{e}^{-y\\ln{z}} = -\\ln{z} $$\r\n%\r\n% $$ -y \\ln{z} = W(-\\ln{z}) $$\r\n%\r\n% $$ y = \\frac{W(-\\ln{z})}{-\\ln{z}} $$\r\n%\r\n% So, the Lambert W function provides a value of the infinite tower\r\n% for any value of $z$.  I didn't know this when I was in junior high school.\r\n\r\n%% The Tower function\r\n% Let's try this formulation for a few values of $z$, including some\r\n% smaller than $\\mbox{e}^{\\mbox{-e}}$.\r\n\r\n\r\n  T = @(z) Lambert_W(-log(z)).\/(-log(z))  \r\n  z = [1.25 sqrt(2) ebar 1.5 e^(-e) 0.05 10.^(-(3:3:15))]';\r\n  disp('         z                    T(z)')\r\n  disp([z T(z)])\r\n\r\n%% Plot\r\n% And now a plot.\r\n\r\n  ezplot(T,[0 ebar])\r\n  axis([0 ebar 0 e])\r\n  title('Tower of Powers')\r\n\r\n%% Unstable\r\n% The value of $y$ obtained from the Lambert W formulation is a fixed point\r\n% of $f(z)$, but a possibly unstable one.  If we have a value of $z$ less\r\n% than $\\mbox{e}^{\\mbox{-e}}$ and we start the fixed point algorithm\r\n% at the value of $y$ given by the Lambert W formula exactly,\r\n% it will stay there.  But if we start with any other value of $y$, the\r\n% iterates will eventually find their way to the pair of bifurcation attractors.  \r\n\r\n\r\n##### SOURCE END ##### e9a2d52dcf9a4d1c856da167477af7f1\r\n-->","protected":false},"excerpt":{"rendered":"<div class=\"overview-image\"><img decoding=\"async\"  class=\"img-responsive\" src=\"https:\/\/blogs.mathworks.com\/images\/cleve\/tower_01.png\" onError=\"this.style.display ='none';\" \/><\/div><!--introduction--><p>I first investigated the tower of powers... <a class=\"read-more\" href=\"https:\/\/blogs.mathworks.com\/cleve\/2013\/09\/16\/tower-of-powers\/\">read more >><\/a><\/p>","protected":false},"author":78,"featured_media":0,"comment_status":"open","ping_status":"closed","sticky":false,"template":"","format":"standard","meta":[],"categories":[12,17],"tags":[],"_links":{"self":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/765"}],"collection":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts"}],"about":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/types\/post"}],"author":[{"embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/users\/78"}],"replies":[{"embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/comments?post=765"}],"version-history":[{"count":5,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/765\/revisions"}],"predecessor-version":[{"id":770,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/765\/revisions\/770"}],"wp:attachment":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/media?parent=765"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/categories?post=765"},{"taxonomy":"post_tag","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/tags?post=765"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}