{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Newton method: Slow convergence" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Consider the function\n", "$$\n", "f(x) = 1/x - a, \\qquad a = 10^{-10}\n", "$$\n", "The root of this function gives the inverse of $a$. We apply Newton method\n", "$$\n", "x_{k+1} = 2 x_k - a x_k^2\n", "$$\n", "starting with initial guess $x_0 = a = 10^{-10}$. First we define the function and its derivative." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "%matplotlib inline\n", "from numpy import abs\n", "\n", "a = 1.0e-10\n", "\n", "def F(x):\n", " return 1/x - a\n", "\n", "def DF(x):\n", " return -1/x**2" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now implement Newton method." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ " 0 2.00000000000000e-10 1.00000000000000e-10 5.00000000000000e+09\n", " 1 4.00000000000000e-10 2.00000000000000e-10 2.50000000000000e+09\n", " 2 8.00000000000000e-10 4.00000000000000e-10 1.25000000000000e+09\n", " 3 1.60000000000000e-09 8.00000000000000e-10 6.25000000000000e+08\n", " 4 3.20000000000000e-09 1.60000000000000e-09 3.12500000000000e+08\n", " 5 6.40000000000000e-09 3.20000000000000e-09 1.56250000000000e+08\n", " 6 1.28000000000000e-08 6.40000000000000e-09 7.81250000000000e+07\n", " 7 2.56000000000000e-08 1.28000000000000e-08 3.90625000000000e+07\n", " 8 5.12000000000000e-08 2.56000000000000e-08 1.95312500000000e+07\n", " 9 1.02400000000000e-07 5.12000000000000e-08 9.76562500000000e+06\n", " 10 2.04800000000000e-07 1.02400000000000e-07 4.88281250000000e+06\n", " 11 4.09600000000000e-07 2.04800000000000e-07 2.44140625000000e+06\n", " 12 8.19200000000000e-07 4.09600000000000e-07 1.22070312500000e+06\n", " 13 1.63840000000000e-06 8.19200000000000e-07 6.10351562500000e+05\n", " 14 3.27680000000000e-06 1.63840000000000e-06 3.05175781250000e+05\n", " 15 6.55360000000000e-06 3.27680000000000e-06 1.52587890625000e+05\n", " 16 1.31072000000000e-05 6.55360000000000e-06 7.62939453124999e+04\n", " 17 2.62144000000000e-05 1.31072000000000e-05 3.81469726562499e+04\n", " 18 5.24287999999999e-05 2.62143999999999e-05 1.90734863281249e+04\n", " 19 1.04857599999999e-04 5.24287999999996e-05 9.53674316406245e+03\n", " 20 2.09715199999998e-04 1.04857599999998e-04 4.76837158203120e+03\n", " 21 4.19430399999991e-04 2.09715199999993e-04 2.38418579101557e+03\n", " 22 8.38860799999965e-04 4.19430399999974e-04 1.19209289550776e+03\n", " 23 1.67772159999986e-03 8.38860799999895e-04 5.96046447753856e+02\n", " 24 3.35544319999944e-03 1.67772159999958e-03 2.98023223876903e+02\n", " 25 6.71088639999775e-03 3.35544319999831e-03 1.49011611938427e+02\n", " 26 1.34217727999910e-02 6.71088639999325e-03 7.45058059691883e+01\n", " 27 2.68435455999640e-02 1.34217727999730e-02 3.72529029845691e+01\n", " 28 5.36870911998559e-02 2.68435455998919e-02 1.86264514922596e+01\n", " 29 1.07374182399424e-01 5.36870911995677e-02 9.31322574610478e+00\n", " 30 2.14748364797694e-01 1.07374182398271e-01 4.65661287302739e+00\n", " 31 4.29496729590777e-01 2.14748364793083e-01 2.32830643648869e+00\n", " 32 8.58993459163107e-01 4.29496729572330e-01 1.16415321821935e+00\n", " 33 1.71798691825243e+00 8.58993459089320e-01 5.82076609084674e-01\n", " 34 3.43597383620971e+00 1.71798691795728e+00 2.91038304517337e-01\n", " 35 6.87194767123882e+00 3.43597383502911e+00 1.45519152233668e-01\n", " 36 1.37438953377553e+01 6.87194766651645e+00 7.27595760918342e-02\n", " 37 2.74877906566211e+01 1.37438953188658e+01 3.63797880209171e-02\n", " 38 5.49755812376843e+01 2.74877905810632e+01 1.81898939854586e-02\n", " 39 1.09951162173137e+02 5.49755809354528e+01 9.09494696772928e-03\n", " 40 2.19902323137349e+02 1.09951160964211e+02 4.54747345886464e-03\n", " 41 4.39804641438994e+02 2.19902318301645e+02 2.27373670443232e-03\n", " 42 8.79609263535176e+02 4.39804622096182e+02 1.13686832721616e-03\n", " 43 1.75921844969911e+03 8.79609186163930e+02 5.68434138608081e-04\n", " 44 3.51843658991326e+03 1.75921814021415e+03 2.84217044304043e-04\n", " 45 7.03687194188691e+03 3.51843535197365e+03 1.42108497152026e-04\n", " 46 1.40737389320171e+04 7.03686699013023e+03 7.10542235760217e-05\n", " 47 2.81474580570215e+04 1.40737191250044e+04 3.55270867880284e-05\n", " 48 5.62948368861036e+04 2.81473788290820e+04 1.77635183940494e-05\n", " 49 1.12589356861341e+05 5.62945199752375e+04 8.88173419709507e-06\n", " 50 2.25177446086354e+05 1.12588089225013e+05 4.44084209868827e-06\n", " 51 4.50349821684486e+05 2.25172375598132e+05 2.22039604962561e-06\n", " 52 9.00679361872783e+05 4.50329540188297e+05 1.11017302537576e-06\n", " 53 1.80127760141428e+06 9.00598239541493e+05 5.55061513813778e-07\n", " 54 3.60223074272882e+06 1.80095314131454e+06 2.77505759158689e-07\n", " 55 7.20316387882525e+06 3.60093313609643e+06 1.38727884082944e-07\n", " 56 1.44011392006640e+07 7.19797532183872e+06 6.93389510486708e-08\n", " 57 2.87815391203003e+07 1.43803999196363e+07 3.46444935387308e-08\n", " 58 5.74802405411872e+07 2.86987014208869e+07 1.72972827981375e-08\n", " 59 1.14630083277107e+08 5.71498427359199e+07 8.62371345646323e-09\n", " 60 2.27946160955002e+08 1.13316077677895e+08 4.28700084182336e-09\n", " 61 4.50696376680593e+08 2.22750215725590e+08 2.11878863851772e-09\n", " 62 8.81080030965884e+08 4.30383654285291e+08 1.03497067786652e-09\n", " 63 1.68452985983508e+09 8.03449828869199e+08 4.93637443801620e-10\n", " 64 3.08529563480257e+09 1.40076577496748e+09 2.24118048435897e-10\n", " 65 5.21868635419195e+09 2.13339071938939e+09 9.16191033777574e-11\n", " 66 7.71390398204098e+09 2.49521762784902e+09 2.96360445149611e-11\n", " 67 9.47737649966719e+09 1.76347251762621e+09 5.51443218860368e-12\n", " 68 9.97268646768999e+09 4.95309968022797e+08 2.73883395397065e-13\n", " 69 9.99992539709528e+09 2.72389294052880e+07 7.46034612872794e-16\n", " 70 9.99999999944344e+09 7.46023481644486e+04 5.56560720338062e-21\n", " 71 1.00000000000000e+10 5.56560720276110e-01 0.00000000000000e+00\n", " 72 1.00000000000000e+10 0.00000000000000e+00 0.00000000000000e+00\n" ] } ], "source": [ "M = 100 # maximum iterations\n", "x = a # initial guess\n", "eps = 1e-15 # relative tolerance on root\n", "\n", "f = F(x)\n", "for i in range(M):\n", " df = DF(x)\n", " dx = -f/df\n", " x = x + dx\n", " e = abs(dx)\n", " f = F(x)\n", " print(\"%6d %22.14e %22.14e %22.14e\" % (i,x,e,abs(f)))\n", " if e < eps * abs(x):\n", " break" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.6.6" } }, "nbformat": 4, "nbformat_minor": 1 }