{"id":2592,"date":"2017-07-17T23:59:48","date_gmt":"2017-07-18T04:59:48","guid":{"rendered":"https:\/\/blogs.mathworks.com\/cleve\/?p=2592"},"modified":"2017-07-18T13:53:57","modified_gmt":"2017-07-18T18:53:57","slug":"what-is-the-condition-number-of-a-matrix","status":"publish","type":"post","link":"https:\/\/blogs.mathworks.com\/cleve\/2017\/07\/17\/what-is-the-condition-number-of-a-matrix\/","title":{"rendered":"What is the Condition Number of a Matrix?"},"content":{"rendered":"<div class=\"content\"><!--introduction--><p>A couple of questions in comments on recent blog posts have prompted me to discuss matrix condition numbers.<\/p><p>In a comment on my post about Hilbert matrices, a reader named Michele asked:<\/p><div><ul><li>Can you comment on when the condition number gives a tight estimate of the error in a computed inverse and whether there is a better estimator?<\/li><\/ul><\/div><p>And in a comment on my post about quadruple precision, Mark asked:<\/p><div><ul><li>Do you have any idea of the slowdown factor ... for large linear equation solves with extremely ill-conditioned matrices?<\/li><\/ul><\/div><p>My short answers are that the error estimate is rarely tight, but that it is not possible to find a better one, and that it takes the same amount of time to solve ill-conditioned linear equations as it does to solve well-conditioned systems.<\/p><!--\/introduction--><h3>Contents<\/h3><div><ul><li><a href=\"#a3219326-029f-4d0b-bf53-12917948c5f2\">Condition Number for Inversion<\/a><\/li><li><a href=\"#0b550be7-f15b-4320-b67a-8f5203e74088\">Norms<\/a><\/li><li><a href=\"#063f47cc-fd5d-40f9-89e0-8963a2c810bb\">Linear equations<\/a><\/li><li><a href=\"#b595d09c-eb46-41b7-ae89-fa58a60cf62e\">Example<\/a><\/li><li><a href=\"#a7ebbbc2-ab18-4fba-b08e-6c72f1cf0ff4\">Close to singular<\/a><\/li><li><a href=\"#b3487b88-ab21-447c-aa7c-bbeeb918a103\">Inverse<\/a><\/li><li><a href=\"#ffac23da-f622-45ba-8107-df5643994e8b\">Speed<\/a><\/li><\/ul><\/div><h4>Condition Number for Inversion<a name=\"a3219326-029f-4d0b-bf53-12917948c5f2\"><\/a><\/h4><p>I should first point out that there are many different condition numbers and that, although the questioners may not have realized it, they were asking about just one of them -- the condition number for inversion.  In general, a condition number applies not only to a particular matrix, but also to the problem being solved.  Are we inverting the matrix, finding its eigenvalues, or computing the exponential?  The list goes on.  A matrix can be poorly conditioned for inversion while the eigenvalue problem is well conditioned. Or, vice versa.<\/p><p>A condition number for a matrix and computational task measures how sensitive the answer is to perturbations in the input data and to roundoff errors made during the solution process.<\/p><p>When we simply say a matrix is \"ill-conditioned\", we are usually just thinking of the sensitivity of its inverse and not of all the other condition numbers.<\/p><h4>Norms<a name=\"0b550be7-f15b-4320-b67a-8f5203e74088\"><\/a><\/h4><p>In order to make these notions more precise, let's start with a vector norm.  Specifically, the <i>Euclidean norm<\/i> or 2- <i>norm<\/i>.<\/p><p>$$\\|x\\| \\ = \\ (\\sum_i x_i^2)^{1\/2}$$<\/p><p>This is the \"as the crow flies\" distance in <i>n<\/i>-dimensional space.<\/p><p>The corresponding norm of a matrix $A$ measures how much the mapping induced by that matrix can stretch vectors.<\/p><p>$$M \\ = \\ \\|A\\| \\ = \\ {\\max {{\\|Ax\\|} \\over {\\|x\\|}}}$$<\/p><p>It is sometimes also important to consider how much a matrix can shrink vectors.<\/p><p>$$m \\ = \\ {\\min {{\\|Ax\\|} \\over {\\|x\\|}}}$$<\/p><p>The reciprocal of the minimum stretching is the norm of the inverse, because<\/p><p>$$m \\ = \\ {\\min {{\\|Ax\\|} \\over {\\|x\\|}}}\r\n\\ = \\ {\\min {{\\|y\\|} \\over {\\|A^{-1} y\\|}}}\r\n\\ = \\ 1\/{\\max {{\\|A^{-1} y\\|} \\over {\\|y\\|}}} \\ = \\ 1\/\\|A^{-1}\\|$$<\/p><p>A <i>singular<\/i> matrix is one that can map nonzero vectors into the zero vector.  For a singular matrix<\/p><p>$$m \\ = \\ 0$$<\/p><p>and the inverse does not exist.<\/p><p>The ratio of the maximum to minimum stretching is the condition number for inversion.<\/p><p>$$\\kappa(A) \\ = \\  {{M} \\over {m}}$$<\/p><p>An equivalent definition is<\/p><p>$$\\kappa(A) \\ = \\  \\|A\\| \\|A^{-1}\\|$$<\/p><p>If a matrix is singular, then its condition number is infinite.<\/p><h4>Linear equations<a name=\"063f47cc-fd5d-40f9-89e0-8963a2c810bb\"><\/a><\/h4><p>The condition number $\\kappa(A)$ is involved in the answer to the question: How much can a change in the right hand side of a system of simultaneous linear equations affect the solution? Consider a system of equations<\/p><p>$$A x \\ = \\ b$$<\/p><p>and a second system obtained by altering the right-hand side<\/p><p>$$A(x + \\delta x) = b + \\delta b$$<\/p><p>Think of $\\delta b$ as being the error in $b$ and $\\delta x$ as being the resulting error in $x$, although we need not make any assumptions that the errors are small.  Because $A (\\delta x) = \\delta b$, the definitions of $M$ and $m$ immediately lead to<\/p><p>$$\\|b\\| \\leq M \\|x\\|$$<\/p><p>and<\/p><p>$$\\|\\delta b\\| \\geq m \\|\\delta x\\|$$<\/p><p>Consequently, if $m \\neq 0$,<\/p><p>$${\\|\\delta x\\| \\over \\|x\\|} \\leq \\kappa(A) {\\|\\delta b\\| \\over\r\n\\|b\\|}$$<\/p><p>The quantity $\\|\\delta b\\|\/\\|b\\|$ is the <i>relative<\/i> change in the right-hand side, and the quantity $\\|\\delta x\\|\/\\|x\\|$ is the resulting <i>relative<\/i> change in the solution.  The advantage of using relative changes is that they are <i>dimensionless<\/i> --- they are not affected by overall scale factors.<\/p><p>This inequality shows that the condition number is a relative error magnification factor. Changes in the right-hand side can cause changes $\\kappa(A)$ times as large in the solution.<\/p><h4>Example<a name=\"b595d09c-eb46-41b7-ae89-fa58a60cf62e\"><\/a><\/h4><p>Let's investigate a system of linear equations involving<\/p><p>$$A=\\left(\\begin{array}{cc} 4.1 &amp; 2.8 \\\\ 9.7 &amp; 6.6 \\end{array}\\right)$$<\/p><p>Take $b$ to be the first column of $A$, so the solution to $Ax = b$ is simply<\/p><p>$$x = \\left(\\begin{array}{c} 1 \\\\ 0 \\end{array}\\right)$$<\/p><p>Switch to MATLAB<\/p><pre class=\"codeinput\">  A = [4.1 2.8; 9.7 6.6]\r\n  b = A(:,1)\r\n  x = A\\b\r\n<\/pre><pre class=\"codeoutput\">A =\r\n    4.1000    2.8000\r\n    9.7000    6.6000\r\nb =\r\n    4.1000\r\n    9.7000\r\nx =\r\n    1.0000\r\n   -0.0000\r\n<\/pre><p>Now add <tt>0.01<\/tt> to the first component of <tt>b<\/tt>.<\/p><pre class=\"codeinput\">  b2 = [4.11; 9.7]\r\n<\/pre><pre class=\"codeoutput\">b2 =\r\n    4.1100\r\n    9.7000\r\n<\/pre><p>The solution changes dramatically.<\/p><pre class=\"codeinput\">  x2 = A\\b2\r\n<\/pre><pre class=\"codeoutput\">x2 =\r\n    0.3400\r\n    0.9700\r\n<\/pre><p>This sensitivity of the solution <tt>x<\/tt> to changes in the right hand side <tt>b<\/tt> is a reflection of the condition number.<\/p><pre class=\"codeinput\">  kappa = cond(A)\r\n<\/pre><pre class=\"codeoutput\">kappa =\r\n   1.6230e+03\r\n<\/pre><p>The upper bound on the possible change in <tt>x<\/tt> indicates changes in all of the significant digits.<\/p><pre class=\"codeinput\">  kappa*norm(b-b2)\/norm(b)\r\n<\/pre><pre class=\"codeoutput\">ans =\r\n    1.5412\r\n<\/pre><p>The actual change in <tt>x<\/tt> resulting from this perturbation is<\/p><pre class=\"codeinput\">  norm(x-x2)\/norm(x)\r\n<\/pre><pre class=\"codeoutput\">ans =\r\n    1.1732\r\n<\/pre><p>So this particular change in the right hand side generated almost the largest possible change in the solution.<\/p><h4>Close to singular<a name=\"a7ebbbc2-ab18-4fba-b08e-6c72f1cf0ff4\"><\/a><\/h4><p>A large condition number means that the matrix is close to being singular.  Let's make a small change in the second row of <tt>A<\/tt>.<\/p><pre class=\"codeinput\">  A\r\n  A2 = [4.1 2.8; 9.676 6.608]\r\n<\/pre><pre class=\"codeoutput\">A =\r\n    4.1000    2.8000\r\n    9.7000    6.6000\r\nA2 =\r\n    4.1000    2.8000\r\n    9.6760    6.6080\r\n<\/pre><p>The resulting matrix is effectively singular.  If we try to compute its inverse, we get a warning message.<\/p><pre class=\"codeinput\">  A2inv = inv(A2)\r\n<\/pre><pre class=\"codeoutput\">Warning: Matrix is close to singular or badly scaled. Results may be inaccurate.\r\nRCOND =  1.988677e-17. \r\nA2inv =\r\n   1.0e+15 *\r\n   -1.4812    0.6276\r\n    2.1689   -0.9190\r\n<\/pre><p>The quantity <tt>RCOND<\/tt> in the warning is an estimate of the reciprocal of the condition number.  Using the reciprocal is left over from the days before we had IEEE floating point arithmetic with <tt>Inf<\/tt> to represent overflow and infinity.  For this example <tt>RCOND<\/tt> is on the order of <tt>eps(1)<\/tt> and the scale factor for <tt>A2inv<\/tt> implies that its elements are useless.  It's impossible to compute something that doesn't exist.<\/p><p>A fairly recent addition to MATLAB is the function <tt>condest<\/tt> that estimates $\\kappa(A)$ .  Now <tt>condest(A)<\/tt> is preferable to <tt>1\/rcond(A)<\/tt>.<\/p><h4>Inverse<a name=\"b3487b88-ab21-447c-aa7c-bbeeb918a103\"><\/a><\/h4><p>The condition number $\\kappa(A)$ also appears in the bound for how much a change $E$ in a matrix $A$ can affect its inverse.<\/p><p>$${{\\|(A+E)^{-1} - A^{-1}\\|} \\over {\\|A^{-1}\\|}} \\ \\le \\\r\n\\kappa(A) {{\\|E\\|} \\over {\\|A\\|}}$$<\/p><p>Jim Wilkinson's work about roundoff error in Gaussian elimination showed that each column of the computed inverse is a column of the exact inverse of a matrix within roundoff error of the given matrix. Let's fudge this a bit and say that <tt>inv(A)<\/tt> computes the exact inverse of $A+E$ where $\\|E\\|$ is on the order of roundoff error compared to $\\|A\\|$ .<\/p><p>We don't know $E$ exactly, but for an <tt>n<\/tt> -by- <tt>n<\/tt> matrix we have the estimate<\/p><p><tt>norm(E)<\/tt> $\\approx \\space$ <tt>n*eps(norm(A))<\/tt><\/p><p>So we have a simple estimate for the error in the computed inverse, relative to the unknown exact inverse.<\/p><p><tt>X =<\/tt> $\\space$ exact inverse of <tt>A<\/tt><\/p><p><tt>Z = inv(A)<\/tt><\/p><p><tt>norm(Z - X)\/norm(X)<\/tt> $\\approx$ <tt>n*eps*cond(A)<\/tt><\/p><p>For our 2-by-2 example the estimate of the relative error in the computed inverse is<\/p><pre class=\"codeinput\">  2*eps*condest(A)\r\n<\/pre><pre class=\"codeoutput\">ans =\r\n   9.9893e-13\r\n<\/pre><p>This says we can expect 12 or 13 (out of 16) significant digits.<\/p><p>Wilkinson had to assume that every individual floating point arithmetic operation incurs the maximum roundoff error.  But only a fraction of the operations have any roundoff error at all and even for those the errors are smaller than the maximum possible.  So this estimate can be expected to an overestimate.  But no tighter estimate is possible.<\/p><p>For our example, the computed inverse is<\/p><pre class=\"codeinput\">  format <span class=\"string\">long<\/span>\r\n  Z = inv(A)\r\n<\/pre><pre class=\"codeoutput\">Z =\r\n -66.000000000000057  28.000000000000025\r\n  97.000000000000085 -41.000000000000036\r\n<\/pre><p>It turns out that the exact inverse has the integer entries produced by<\/p><pre class=\"codeinput\">  X = round(Z)\r\n<\/pre><pre class=\"codeoutput\">X =\r\n   -66    28\r\n    97   -41\r\n<\/pre><p>We can compare the actual relative error with the estimate.<\/p><pre class=\"codeinput\">  format <span class=\"string\">short<\/span>\r\n  norm(Z - X)\/norm(X)\r\n<\/pre><pre class=\"codeoutput\">ans =\r\n   8.7342e-16\r\n<\/pre><p>So we actually have more than 15 significant digits of accuracy.<\/p><h4>Speed<a name=\"ffac23da-f622-45ba-8107-df5643994e8b\"><\/a><\/h4><p>In simplified outline, the algorithm for computing the inverse of an $n$ -by- $n$ matrix, or for solving a system of $n$ linear equations, involves loops of length $n$ nested three deep.  Each of the $n^2$ elements is accessed roughly $n$ times.  So the computational complexity is proportional to $n^3$.  Lots of things can be done that affect the coefficient of proportionality, but it's still order $n^3$.<\/p><p>The actual numbers in the matrix (generally) don't affect the execution time.  A nearly singular matrix can be inverted just as fast as a well-conditioned one.  The answer might not be very accurate if the condition number is large, but $\\kappa(A)$ does not play a role in the speed.<\/p><script language=\"JavaScript\"> <!-- \r\n    function grabCode_db324eef4c604f70bfe051afe8b3fd25() {\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='db324eef4c604f70bfe051afe8b3fd25 ' + '##### ' + 'SOURCE BEGIN' + ' #####';\r\n        t2='##### ' + 'SOURCE END' + ' #####' + ' db324eef4c604f70bfe051afe8b3fd25';\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 2017 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_db324eef4c604f70bfe051afe8b3fd25()\"><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; R2017a<br><\/p><\/div><!--\r\ndb324eef4c604f70bfe051afe8b3fd25 ##### SOURCE BEGIN #####\r\n%% What is the Condition Number of a Matrix?\r\n% A couple of questions in comments on recent blog posts\r\n% have prompted me to discuss matrix condition numbers.\r\n%\r\n% In a comment on my post about Hilbert matrices, a reader\r\n% named Michele asked:\r\n%\r\n% * Can you comment on when the condition number gives a tight estimate \r\n% of the error in a computed inverse and whether there is a better\r\n% estimator?\r\n%\r\n% And in a comment on my post about quadruple precision, Mark asked:\r\n%\r\n% * Do you have any idea of the slowdown factor ... for large\r\n% linear equation solves with extremely ill-conditioned matrices?\r\n%\r\n% My short answers are that the error estimate is rarely tight, but that\r\n% it is not possible to find a better one, and that it takes the same\r\n% amount of time to solve ill-conditioned linear equations as it does\r\n% to solve well-conditioned systems.\r\n\r\n%% Condition Number for Inversion\r\n% I should first point out that there are many different condition numbers\r\n% and that, although the questioners may not have realized it, they\r\n% were asking about just one of them REPLACE_WITH_DASH_DASH the condition number for\r\n% inversion.  In general, a condition number applies not only to a\r\n% particular matrix, but also to the problem being solved.  Are we\r\n% inverting the matrix, finding its eigenvalues, or computing the\r\n% exponential?  The list goes on.  A matrix can be poorly conditioned\r\n% for inversion while the eigenvalue problem is well conditioned.\r\n% Or, vice versa.\r\n\r\n%%\r\n% A condition number for a matrix and computational task measures how\r\n% sensitive the answer is to perturbations in the input data and to\r\n% roundoff errors made during the solution process.\r\n\r\n%%\r\n% When we simply say a matrix is \"ill-conditioned\", we are usually just\r\n% thinking of the sensitivity of its inverse and not of all the other\r\n% condition numbers.\r\n\r\n%% Norms\r\n% In order to make these notions more precise, let's start with a\r\n% vector norm.  Specifically, the _Euclidean norm_ or 2- _norm_.\r\n%\r\n% $$\\|x\\| \\ = \\ (\\sum_i x_i^2)^{1\/2}$$\r\n%\r\n% This is the \"as the crow flies\" distance in _n_-dimensional space.\r\n\r\n%%\r\n% The corresponding norm of a matrix $A$ measures how much the mapping\r\n% induced by that matrix can stretch vectors.\r\n%\r\n% $$M \\ = \\ \\|A\\| \\ = \\ {\\max {{\\|Ax\\|} \\over {\\|x\\|}}}$$\r\n%\r\n% It is sometimes also important to consider how much a matrix can shrink\r\n% vectors.\r\n%\r\n% $$m \\ = \\ {\\min {{\\|Ax\\|} \\over {\\|x\\|}}}$$\r\n\r\n%%\r\n% The reciprocal of the minimum stretching\r\n% is the norm of the inverse, because\r\n%\r\n% $$m \\ = \\ {\\min {{\\|Ax\\|} \\over {\\|x\\|}}}\r\n% \\ = \\ {\\min {{\\|y\\|} \\over {\\|A^{-1} y\\|}}}   \r\n% \\ = \\ 1\/{\\max {{\\|A^{-1} y\\|} \\over {\\|y\\|}}} \\ = \\ 1\/\\|A^{-1}\\|$$\r\n\r\n%%\r\n% A _singular_ matrix is one that can map nonzero vectors into the\r\n% zero vector.  For a singular matrix\r\n%\r\n% $$m \\ = \\ 0$$\r\n%\r\n% and the inverse does not exist.\r\n%%\r\n% The ratio of the maximum to minimum stretching is the condition number\r\n% for inversion.\r\n%\r\n% $$\\kappa(A) \\ = \\  {{M} \\over {m}}$$\r\n\r\n%%\r\n% An equivalent definition is\r\n%\r\n% $$\\kappa(A) \\ = \\  \\|A\\| \\|A^{-1}\\|$$\r\n\r\n%%\r\n% If a matrix is singular, then its condition number is infinite.\r\n\r\n%% Linear equations\r\n% The condition number $\\kappa(A)$ is involved in the answer\r\n% to the question: \r\n% How much can a change in the right hand side of a system\r\n% of simultaneous linear equations affect the solution?\r\n% Consider a system of equations\r\n%\r\n% $$A x \\ = \\ b$$\r\n%\r\n% and a second system obtained by altering the right-hand side\r\n%\r\n% $$A(x + \\delta x) = b + \\delta b$$\r\n%\r\n% Think of $\\delta b$ as being the error in $b$ and $\\delta x$ as being\r\n% the resulting error in $x$, although we need not make any assumptions \r\n% that the errors are small.  Because $A (\\delta x) = \\delta b$, the\r\n% definitions of $M$ and $m$ immediately lead to\r\n%\r\n% $$\\|b\\| \\leq M \\|x\\|$$\r\n% \r\n% and\r\n%\r\n% $$\\|\\delta b\\| \\geq m \\|\\delta x\\|$$\r\n% \r\n% Consequently, if $m \\neq 0$,\r\n% \r\n% $${\\|\\delta x\\| \\over \\|x\\|} \\leq \\kappa(A) {\\|\\delta b\\| \\over \r\n% \\|b\\|}$$\r\n\r\n%%\r\n% The quantity $\\|\\delta b\\|\/\\|b\\|$ is the _relative_ change in the\r\n% right-hand side, and the quantity $\\|\\delta x\\|\/\\|x\\|$ is the resulting\r\n% _relative_ change in the solution.  The advantage of using relative\r\n% changes is that they are _dimensionless_ REPLACE_WITH_DASH_DASH- they are not affected by\r\n% overall scale factors.\r\n\r\n%%\r\n% This inequality shows that the condition number is a relative error \r\n% magnification factor.\r\n% Changes in the right-hand side can cause changes $\\kappa(A)$ times as \r\n% large in the solution.\r\n\r\n%% Example\r\n% Let's investigate a system of linear equations involving\r\n%\r\n% $$A=\\left(\\begin{array}{cc} 4.1 & 2.8 \\\\ 9.7 & 6.6 \\end{array}\\right)$$\r\n%\r\n% Take $b$ to be the first column of $A$, so the solution to $Ax = b$\r\n% is simply\r\n%\r\n% $$x = \\left(\\begin{array}{c} 1 \\\\ 0 \\end{array}\\right)$$\r\n\r\n%%\r\n% Switch to MATLAB\r\n\r\n  A = [4.1 2.8; 9.7 6.6]\r\n  b = A(:,1)\r\n  x = A\\b\r\n  \r\n%%\r\n% Now add |0.01| to the first component of |b|.\r\n%\r\n  b2 = [4.11; 9.7]\r\n  \r\n%%\r\n% The solution changes dramatically.\r\n\r\n  x2 = A\\b2\r\n  \r\n%%\r\n% This sensitivity of the solution |x| to changes in the right hand side\r\n% |b| is a reflection of the condition number.\r\n\r\n  kappa = cond(A)\r\n  \r\n%%\r\n% The upper bound on the possible change in |x| indicates changes\r\n% in all of the significant digits.\r\n\r\n  kappa*norm(b-b2)\/norm(b)\r\n  \r\n%%\r\n% The actual change in |x| resulting from this perturbation is\r\n\r\n  norm(x-x2)\/norm(x)\r\n  \r\n%%\r\n% So this particular change in the right hand side generated almost the\r\n% largest possible change in the solution.\r\n   \r\n%% Close to singular\r\n% A large condition number means that the matrix is close to being\r\n% singular.  Let's make a small change in the second row of |A|.\r\n\r\n  A\r\n  A2 = [4.1 2.8; 9.676 6.608]\r\n  \r\n%%\r\n% The resulting matrix is effectively singular.  If we try to compute\r\n% its inverse, we get a warning message.\r\n\r\n  A2inv = inv(A2)\r\n  \r\n%%\r\n% The quantity |RCOND| in the warning is an estimate of the reciprocal of \r\n% the condition number.  Using the reciprocal is left over from the days \r\n% before we had IEEE floating point arithmetic with |Inf| to represent \r\n% overflow and infinity.  For this example |RCOND| is on the order of\r\n% |eps(1)| and the scale factor for |A2inv| implies that its elements \r\n% are useless.  It's impossible to compute something that doesn't exist.\r\n\r\n%%\r\n% A fairly recent addition to MATLAB is the function |condest| that\r\n% estimates $\\kappa(A)$ .  Now |condest(A)| is preferable to |1\/rcond(A)|.\r\n\r\n%% Inverse\r\n% The condition number $\\kappa(A)$ also appears in the bound for how\r\n% much a change $E$ in a matrix $A$ can affect its inverse.\r\n%\r\n% $${{\\|(A+E)^{-1} - A^{-1}\\|} \\over {\\|A^{-1}\\|}} \\ \\le \\\r\n% \\kappa(A) {{\\|E\\|} \\over {\\|A\\|}}$$\r\n\r\n%%\r\n% Jim Wilkinson's work about roundoff error in Gaussian elimination\r\n% showed that each column of the computed inverse is a column of\r\n% the exact inverse of a matrix within roundoff error of the given matrix.\r\n% Let's fudge this a bit and say that |inv(A)| computes the exact inverse\r\n% of $A+E$ where $\\|E\\|$ is on the order of roundoff error compared to\r\n% $\\|A\\|$ .\r\n\r\n%%\r\n% We don't know $E$ exactly, but for an |n| -by- |n| matrix\r\n% we have the estimate\r\n%\r\n% |norm(E)| $\\approx \\space$ |n*eps(norm(A))|\r\n%\r\n% So we have a simple estimate for the error in the computed inverse,\r\n% relative to the unknown exact inverse.\r\n%\r\n% |X =| $\\space$ exact inverse of |A|\r\n%\r\n% |Z = inv(A)|\r\n%\r\n% |norm(Z - X)\/norm(X)| $\\approx$ |n*eps*cond(A)|\r\n\r\n%%\r\n% For our 2-by-2 example the estimate of the relative error in the\r\n% computed inverse is\r\n\r\n  2*eps*condest(A)\r\n\r\n%%\r\n% This says we can expect 12 or 13 (out of 16) significant digits.\r\n\r\n%%\r\n% Wilkinson had to assume that every individual floating point arithmetic\r\n% operation incurs the maximum roundoff error.  But only a fraction of the\r\n% operations have any roundoff error at all and even for those the errors\r\n% are smaller than the maximum possible.  So this estimate can be\r\n% expected to an overestimate.  But no tighter estimate is possible.\r\n\r\n%%\r\n% For our example, the computed inverse is\r\n\r\n  format long\r\n  Z = inv(A)\r\n  \r\n%%\r\n% It turns out that the exact inverse has the integer entries produced by\r\n\r\n  X = round(Z)\r\n  \r\n%%\r\n% We can compare the actual relative error with the estimate.\r\n\r\n  format short\r\n  norm(Z - X)\/norm(X)\r\n  \r\n%%\r\n% So we actually have more than 15 significant digits of accuracy.\r\n\r\n%% Speed\r\n% In simplified outline, the algorithm for computing the inverse of an\r\n% $n$ -by- $n$ matrix, or for solving a system of $n$ linear equations,\r\n% involves loops of length $n$ nested three deep.  Each of the $n^2$\r\n% elements is accessed roughly $n$ times.  So the computational complexity\r\n% is proportional to $n^3$.  Lots of things can be done that affect the\r\n% coefficient of proportionality, but it's still order $n^3$.\r\n\r\n%%\r\n% The actual numbers in the matrix (generally) don't affect the execution\r\n% time.  A nearly singular matrix can be inverted just as fast as a \r\n% well-conditioned one.  The answer might not be very accurate if\r\n% the condition number is large, but $\\kappa(A)$ does not play a role\r\n% in the speed.\r\n##### SOURCE END ##### db324eef4c604f70bfe051afe8b3fd25\r\n-->","protected":false},"excerpt":{"rendered":"<div class=\"overview-image\"><img src=\"https:\/\/blogs.mathworks.com\/cleve\/files\/condition_thumb-1.png\" class=\"img-responsive attachment-post-thumbnail size-post-thumbnail wp-post-image\" alt=\"\" decoding=\"async\" loading=\"lazy\" \/><\/div><!--introduction--><p>A couple of questions in comments on recent blog posts have prompted me to discuss matrix condition numbers.... <a class=\"read-more\" href=\"https:\/\/blogs.mathworks.com\/cleve\/2017\/07\/17\/what-is-the-condition-number-of-a-matrix\/\">read more >><\/a><\/p>","protected":false},"author":78,"featured_media":2620,"comment_status":"open","ping_status":"closed","sticky":false,"template":"","format":"standard","meta":[],"categories":[6,16,14],"tags":[],"_links":{"self":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/2592"}],"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=2592"}],"version-history":[{"count":1,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/2592\/revisions"}],"predecessor-version":[{"id":2593,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/2592\/revisions\/2593"}],"wp:featuredmedia":[{"embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/media\/2620"}],"wp:attachment":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/media?parent=2592"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/categories?post=2592"},{"taxonomy":"post_tag","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/tags?post=2592"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}