{"id":385,"date":"2012-11-12T12:30:23","date_gmt":"2012-11-12T17:30:23","guid":{"rendered":"https:\/\/blogs.mathworks.com\/cleve\/?p=385"},"modified":"2016-12-05T13:55:30","modified_gmt":"2016-12-05T18:55:30","slug":"magic-squares-part-3-linear-algebra","status":"publish","type":"post","link":"https:\/\/blogs.mathworks.com\/cleve\/2012\/11\/12\/magic-squares-part-3-linear-algebra\/","title":{"rendered":"Magic Squares, Part 3, Linear Algebra"},"content":{"rendered":"&nbsp;\r\n<div class=\"content\"><!--introduction-->We know a little bit, but not much, about the linear algebraic properties of magic squares.\r\n\r\n<!--\/introduction-->\r\n<h3>Contents<\/h3>\r\n<div>\r\n<ul>\r\n \t<li><a href=\"#444a61af-8bac-42a1-9401-04bb846df002\">Rank<\/a><\/li>\r\n \t<li><a href=\"#6c162ac2-9ddc-43b2-b308-318aea53da64\">Odd Order<\/a><\/li>\r\n \t<li><a href=\"#b4f2bfcb-31e9-42d9-889f-c4572d7abfff\">Doubly Even Order<\/a><\/li>\r\n \t<li><a href=\"#067ff3a2-a52b-4707-9331-afa766551d28\">Order Four<\/a><\/li>\r\n \t<li><a href=\"#f93daa68-6669-439a-b510-4521d07ac5da\">Singly Even<\/a><\/li>\r\n \t<li><a href=\"#794eca4b-8b9c-4d2b-9005-702d3fe2a063\">Further Reading<\/a><\/li>\r\n<\/ul>\r\n<\/div>\r\n<h4>Rank<a name=\"444a61af-8bac-42a1-9401-04bb846df002\"><\/a><\/h4>\r\nHere is my favorite experiment involving the linear algebraic properties of magic squares. Recall that the <i>rank<\/i> of a matrix is the number of linearly independent rows or columns. If an <i>n<\/i> -by- <i>n<\/i> matrix has <i>full rank<\/i>, that is rank equal to <i>n<\/i>, then it is <i>nonsingular<\/i> and has an inverse. Let's look at the rank of the magic squares generated by MATLAB.\r\n<pre class=\"codeinput\">    r = zeros(1,33);\r\n    <span class=\"keyword\">for<\/span> n = 3:33\r\n       r(n) = rank(magic(n));\r\n    <span class=\"keyword\">end<\/span>\r\n    bar(r)\r\n    xlabel(<span class=\"string\">'n'<\/span>)\r\n    ylabel(<span class=\"string\">'rank'<\/span>)\r\n    title(<span class=\"string\">'rank(magic(n))'<\/span>)\r\n    snapnow\r\n<\/pre>\r\n<img decoding=\"async\" src=\"https:\/\/blogs.mathworks.com\/images\/cleve\/magic_3_01.png\" alt=\"\" hspace=\"5\" vspace=\"5\" \/>\r\n\r\nWe can clearly see three different cases. Odd order magic squares have full rank, singly even order magic squares have rank equal to roughly half their order, and doubly even order magic squares of any order have rank equal to three. Why is that? It reflects the three different cases in the MATLAB magic square generator.\r\n<h4>Odd Order<a name=\"6c162ac2-9ddc-43b2-b308-318aea53da64\"><\/a><\/h4>\r\nRecall from last week that we have two different ways of generating the odd order magic square produced by <tt>magic(n)<\/tt>. De La Loubere's method generates the matrix elements one-by-one along broken diagonals in a manner that seems very unlikely to generate linear dependence between rows or columns.\r\n\r\nOur new algorithm generates the same magic square with a few array operations.\r\n<pre>   [I,J] = ndgrid(1:n);\r\n   A = mod(I+J+(n-3)\/2,n);\r\n   B = mod(I+2*J-2,n);\r\n   M = n*A + B + 1;<\/pre>\r\nIt is not hard to prove that <tt>A<\/tt> and <tt>B<\/tt> have full rank and so it is likely that their sum has full rank as well.\r\n\r\nThere are many other magic squares besides the ones generated by MATLAB. Are all the odd ordered ones nonsingular? Frankly, I have no idea. The magic square property tells us very little about linear independence.\r\n\r\nDoes anybody know of a 5-by-5 magic square that is also a singular matrix? I would love to see it. Or, can anybody prove that odd-ordered magic squares must be nonsingular? I would love to see a proof as well. Actually, I better see just one or the other.\r\n<h4>Doubly Even Order<a name=\"b4f2bfcb-31e9-42d9-889f-c4572d7abfff\"><\/a><\/h4>\r\nFor doubly even order, the rank is always three. When I first saw that -- years ago -- I was really surprised. For a long time I had a vague idea that the rank could be explained by the simple algorithm that generates these matrices. But now I'm not so sure.\r\n\r\nLet's look at the generation of <tt>magic(4)<\/tt>. Start with the integers <tt>1:16<\/tt> arranged in a 4-by-4 array. This, of course, is a matrix with rank 1, but we're about to alter it drastically. Mark half of the elements with an asterisk in a pattern like this.\r\n<pre>   1*    2     3     4*\r\n   5     6*    7*    8\r\n   9    10*   11*   12\r\n  13*   14    15    16*<\/pre>\r\nNow regard the marked elements as a one-dimensional array and flip that array end-for-end. This produces <tt>magic(4)<\/tt>.\r\n<pre>  16*    2     3    13*\r\n   5    11*   10*    8\r\n   9     7*    6*   12\r\n   4*   14    15     1*<\/pre>\r\nIt must be easy to describe that indexing operation in matrix terms in a way that generalizes to higher order and explains why the resulting matrix happens to have rank 3. But I don't see how today.\r\n<h4>Order Four<a name=\"067ff3a2-a52b-4707-9331-afa766551d28\"><\/a><\/h4>\r\nWhat about other 4-by-4 magic squares? Do they all have rank 3? Are there any nonsingular magic squares of order 4? I am able to answer that question because I came across this terrific paper by Bill Press, one of the authors of <i>Numerical Recipes<\/i>. <a href=\"\">Did D\u00fcrer Intentionally Show Only His Second-Best Magic Square?<\/a> (I won't spoil Bill's punch line by revealing the answer to the question he poises in his title.)\r\n\r\nIn order to investigate D\u00fcrer's possible choices for a special magic square, Press describes a seven variable parameterization of the 16 elements in a 4-by-4 array that characterizes all possible magic squares of order 4. (He attributes the parameterization to Maurice Kraitchik.) This approach leads to a program, <a href=\"https:\/\/blogs.mathworks.com\/images\/cleve\/enumerate4.m\">enumerate.m<\/a>, that checks the rank of all magic squares of order 4. This is the first time I have ever written a MATLAB code with <tt>for<\/tt> loops nested seven deep. It generates\r\n<pre class=\"language-matlab\">16*15*14*13*12*11*10\r\n<\/pre>\r\n<pre class=\"language-matlab\">= 57657600\r\n<\/pre>\r\nmatrices. Of these, only 7040 pass the check in the inner block to qualify as legal magic squares.\r\n\r\nDividing the final counts by 8 to account for eight-fold symmetry, we find that that there are 880 magic squares of order four. Their ranks are\r\n<pre>  rank  count\r\n    3    640\r\n    4    240<\/pre>\r\nSo, there ware nonsingular 4-by-4 magic squares. Other than that, I attach no great significant to the actual counts. The process of finding them is more interesting than the counts themselves.\r\n\r\nBy the way, the entire computation by <tt>enumerate4<\/tt> takes about 415 seconds on my Lenovo X201 2.67 GHz i7 laptop. The first nonsingular 4-by-4 magic square that <tt>enumerate4<\/tt> finds is\r\n<pre>   1     2    16    15\r\n  13    14     4     3\r\n  12     7     9     6\r\n   8    11     5    10<\/pre>\r\n<h4>Singly Even<a name=\"f93daa68-6669-439a-b510-4521d07ac5da\"><\/a><\/h4>\r\nThe bar chart above shows that the rank for singly even magic squares, where <i>n<\/i> is divisible by two but not by four, is about $n\/2$. Actually, the rank is $n\/2+2$, and it is not hard to see why.\r\n\r\nThe first step of the singly even algorithm is to stack together four copies of the square of odd order $n\/2$.\r\n<pre>  p = n\/2;\r\n  A = magic(p);\r\n  M = [A A+2*p^2; A+3*p^2 A+p^2];<\/pre>\r\nThis <tt>M<\/tt> has equal column sums, but not equal row sums. And, since <tt>A<\/tt> has rank $n\/2$, this <tt>M<\/tt> has rank $n\/2+1$\r\n\r\nThe next step interchanges some blocks of elements to correct the row sums.\r\n<pre>  i = (1:p)';\r\n  k = (n-2)\/4;\r\n  j = [1:k (n-k+2):n];\r\n  M([i; i+p],j) = M([i+p; i],j);<\/pre>\r\nIt turns out that this step does not change the rank. Now this <tt>M<\/tt> has equal row and column sums, but its diagonal sums are not quite right.\r\n\r\nThe final step swaps two pairs of elements to correct the diagonal.\r\n<pre>  i = k+1;\r\n  j = [1 i];\r\n  M([i; i+p],j) = M([i+p; i],j);<\/pre>\r\nThis is a rank one change that produces the final magic square and increases the rank to $n\/2+2$.\r\n<h4>Further Reading<a name=\"794eca4b-8b9c-4d2b-9005-702d3fe2a063\"><\/a><\/h4>\r\nWilliam H. Press, <a href=\"\">Did D\u00fcrer Intentionally Show Only His Second-Best Magic Square?<\/a>, 2009.\r\n\r\nCBM, <a href=\"https:\/\/www.mathworks.com\/content\/dam\/mathworks\/mathworks-dot-com\/moler\/exm\/chapters\/magic.pdf\">Magic Squares<\/a> chapter of <a href=\"https:\/\/www.mathworks.com\/moler\/exm\"><i>Experiments with MATLAB<\/i><\/a>, 2011.\r\n\r\n<script>\/\/ <![CDATA[\r\nfunction grabCode_5681b7a27e4a4cbe81cadb7ca43801c0() {\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='5681b7a27e4a4cbe81cadb7ca43801c0 ' + '##### ' + 'SOURCE BEGIN' + ' #####';\r\n        t2='##### ' + 'SOURCE END' + ' #####' + ' 5681b7a27e4a4cbe81cadb7ca43801c0';\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 2012 The MathWorks, Inc.';\r\n\r\n        w = window.open();\r\n        d = w.document;\r\n        d.write('\r\n\r\n<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>\r\n\r\n\r\n\\n');\r\n\r\n        d.title = title + ' (MATLAB code)';\r\n        d.close();\r\n    }\r\n\/\/ ]]><\/script>\r\n<p style=\"text-align: right; font-size: xx-small; font-weight: lighter; font-style: italic; color: gray;\">\r\n<a><span style=\"font-size: x-small; font-style: italic;\">Get\r\nthe MATLAB code<noscript>(requires JavaScript)<\/noscript><\/span><\/a>\r\n\r\nPublished with MATLAB\u00ae R2012b<\/p>\r\n<p class=\"footer\">\r\nPublished with MATLAB\u00ae R2012b<\/p>\r\n\r\n<\/div>\r\n<!--\r\n5681b7a27e4a4cbe81cadb7ca43801c0 ##### SOURCE BEGIN #####\r\n%% Magic Squares, Part 3, Linear Algebra\r\n% We know a little bit, but not much, about the linear algebraic\r\n% properties of magic squares.\r\n\r\n%% Rank\r\n% Here is my favorite experiment involving the linear algebraic\r\n% properties of magic squares.  Recall that the _rank_ of a matrix\r\n% is the number of linearly independent rows or columns.\r\n% If an _n_ -by- _n_ matrix has _full rank_, that is rank equal to _n_,\r\n% then it is _nonsingular_ and has an inverse.  Let's look at the rank\r\n% of the magic squares generated by MATLAB.\r\n\r\nr = zeros(1,33);\r\nfor n = 3:33\r\nr(n) = rank(magic(n));\r\nend\r\nbar(r)\r\nxlabel('n')\r\nylabel('rank')\r\ntitle('rank(magic(n))')\r\nsnapnow\r\n\r\n%%\r\n% We can clearly see three different cases.  Odd order magic squares\r\n% have full rank, singly even order magic squares have rank equal to\r\n% roughly half their order, and doubly even order magic squares of any order\r\n% have rank equal to three.  Why is that?  It reflects the three different\r\n% cases in the MATLAB magic square generator.\r\n\r\n%% Odd Order\r\n% Recall from last week that we have two different ways of generating\r\n% the odd order magic square produced by |magic(n)|.\r\n% De La Loubere's method generates the matrix elements one-by-one along\r\n% broken diagonals in a manner that seems very unlikely to generate\r\n% linear dependence between rows or columns.\r\n%\r\n% Our new algorithm generates the same magic square with a few array\r\n% operations.\r\n%\r\n%     [I,J] = ndgrid(1:n);\r\n%     A = mod(I+J+(n-3)\/2,n);\r\n%     B = mod(I+2*J-2,n);\r\n%     M = n*A + B + 1;\r\n%\r\n% It is not hard to prove that |A| and |B| have full rank and so it is\r\n% likely that their sum has full rank as well.\r\n%\r\n% There are many other magic squares besides the ones generated by MATLAB.\r\n% Are all the odd ordered ones nonsingular?  Frankly, I have no idea.\r\n% The magic square property tells us very little about linear independence.\r\n%\r\n% Does anybody know of a 5-by-5 magic square that is also a singular matrix?\r\n% I would love to see it.\r\n% Or, can anybody prove that odd-ordered magic squares must be nonsingular?\r\n% I would love to see a proof as well.\r\n% Actually, I better see just one or the other.\r\n\r\n%% Doubly Even Order\r\n% For doubly even order, the rank is always three.  When I first saw that\r\n% REPLACE_WITH_DASH_DASH years ago REPLACE_WITH_DASH_DASH I was really surprised.  For a long time I had a vague idea\r\n% that the rank could be explained by the simple algorithm that generates\r\n% these matrices.  But now I'm not so sure.\r\n%\r\n% Let's look at the generation of |magic(4)|.\r\n% Start with the integers |1:16| arranged in a 4-by-4 array.\r\n% This, of course, is a matrix with rank 1, but we're about to alter it\r\n% drastically.\r\n% Mark half of the elements with an asterisk in a pattern like this.\r\n%\r\n%     1*    2     3     4*\r\n%     5     6*    7*    8\r\n%     9    10*   11*   12\r\n%    13*   14    15    16*\r\n%\r\n% Now regard the marked elements as a one-dimensional array and\r\n% flip that array end-for-end.  This produces |magic(4)|.\r\n%\r\n%    16*    2     3    13*\r\n%     5    11*   10*    8\r\n%     9     7*    6*   12\r\n%     4*   14    15     1*\r\n%\r\n% It must be easy to describe that indexing operation in matrix terms\r\n% in a way that generalizes to higher order and explains\r\n% why the resulting matrix happens to have rank 3.\r\n% But I don't see how today.\r\n\r\n%% Order Four\r\n% What about other 4-by-4 magic squares? Do they all have rank 3?\r\n% Are there any nonsingular magic squares of order 4?\r\n% I am able to answer that question because\r\n% I came across this terrific paper by Bill Press, one of the authors of\r\n% _Numerical Recipes_.  < % Did D\u251c\u255drer Intentionally Show Only His Second-Best Magic Square?>\r\n% (I won't spoil Bill's punch line by revealing the answer to the question\r\n% he poises in his title.)\r\n%\r\n% In order to investigate D\u251c\u255drer's possible choices for a special magic square,\r\n% Press describes a seven variable parameterization of the 16 elements\r\n% in a 4-by-4 array that characterizes all possible magic squares of order 4.\r\n% (He attributes the parameterization to Maurice Kraitchik.)\r\n% This approach leads to a program,\r\n% <https:\/\/blogs.mathworks.com\/images\/cleve\/enumerate4.m enumerate.m>,\r\n% that checks the rank of all magic squares of order 4.\r\n% This is the first time I have ever written a MATLAB code\r\n% with |for| loops nested seven deep.  It generates\r\n%\r\n%   16*15*14*13*12*11*10\r\n%\r\n%   = 57657600\r\n%\r\n% matrices.  Of these, only 7040 pass the check in the inner block to qualify\r\n% as legal magic squares.\r\n%\r\n% Dividing the final counts by 8 to account for eight-fold symmetry,\r\n% we find that that there are 880 magic squares of order four.\r\n% Their ranks are\r\n%\r\n%    rank  count\r\n%      3    640\r\n%      4    240\r\n%\r\n% So, there ware nonsingular 4-by-4 magic squares.\r\n% Other than that, I attach no great significant to the actual counts.\r\n% The process of finding them is more interesting than the counts themselves.\r\n%\r\n% By the way, the entire computation by |enumerate4| takes about 415 seconds\r\n% on my Lenovo X201 2.67 GHz i7 laptop.\r\n% The first nonsingular 4-by-4 magic square that |enumerate4| finds is\r\n%\r\n%     1     2    16    15\r\n%    13    14     4     3\r\n%    12     7     9     6\r\n%     8    11     5    10\r\n\r\n%% Singly Even\r\n% The bar chart above shows that the rank for singly even magic squares,\r\n% where _n_ is divisible by two but not by four, is about $n\/2$.\r\n% Actually, the rank is $n\/2+2$, and it is not hard to see why.\r\n%\r\n% The first step of the singly even algorithm is to stack together\r\n% four copies of the square of odd order $n\/2$.\r\n%\r\n%    p = n\/2;\r\n%    A = magic(p);\r\n%    M = [A A+2*p^2; A+3*p^2 A+p^2];\r\n%\r\n% This |M| has equal column sums, but not equal row sums.  And, since |A| has\r\n% rank $n\/2$, this |M| has rank $n\/2+1$\r\n%\r\n% The next step interchanges some blocks of elements to correct the row sums.\r\n%\r\n%    i = (1:p)';\r\n%    k = (n-2)\/4;\r\n%    j = [1:k (n-k+2):n];\r\n%    M([i; i+p],j) = M([i+p; i],j);\r\n%\r\n% It turns out that this step does not change the rank.  Now this |M| has\r\n% equal row and column sums, but its diagonal sums are not quite right.\r\n%\r\n% The final step swaps two pairs of elements to correct the diagonal.\r\n%\r\n%    i = k+1;\r\n%    j = [1 i];\r\n%    M([i; i+p],j) = M([i+p; i],j);\r\n%\r\n% This is a rank one change that produces the final magic square and\r\n% increases the rank to $n\/2+2$.\r\n\r\n%% Further Reading\r\n%\r\n% William H. Press, < % Did D\u251c\u255drer Intentionally Show Only His Second-Best Magic Square?>,\r\n% 2009.\r\n%\r\n% CBM,\r\n% <https:\/\/www.mathworks.com\/moler\/exm\/chapters\/magic.pdf Magic Squares>\r\n% chapter of\r\n% <https:\/\/www.mathworks.com\/moler\/exm _Experiments with MATLAB_>, 2011.\r\n\r\n##### SOURCE END ##### 5681b7a27e4a4cbe81cadb7ca43801c0\r\n-->","protected":false},"excerpt":{"rendered":"<!--introduction-->We know a little bit, but not much, about the linear algebraic properties of magic squares.\r\n\r\n<!--\/introduction-->... <a class=\"read-more\" href=\"https:\/\/blogs.mathworks.com\/cleve\/2012\/11\/12\/magic-squares-part-3-linear-algebra\/\">read more >><\/a><\/p>","protected":false},"author":78,"featured_media":0,"comment_status":"open","ping_status":"closed","sticky":false,"template":"","format":"standard","meta":[],"categories":[9],"tags":[],"_links":{"self":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/385"}],"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=385"}],"version-history":[{"count":9,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/385\/revisions"}],"predecessor-version":[{"id":2179,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/385\/revisions\/2179"}],"wp:attachment":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/media?parent=385"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/categories?post=385"},{"taxonomy":"post_tag","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/tags?post=385"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}