Cleve Moler is the author of the first MATLAB, one of the founders of MathWorks, and is currently Chief Mathematician at the company. He is the author of two books about MATLAB that are available online. He writes here about MATLAB, scientific computing and interesting mathematics.
Here is my favorite experiment involving the linear algebraic properties of magic squares. Recall that the rank of a matrix is the number of linearly independent rows or columns. If an n -by- n matrix has full rank, that is rank equal to n, then it is nonsingular and has an inverse. Let's look at the rank of the magic squares generated by MATLAB.
r = zeros(1,33);
for n = 3:33
r(n) = rank(magic(n));
end
bar(r)
xlabel('n')
ylabel('rank')
title('rank(magic(n))')
snapnow
We 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.
Odd Order
Recall from last week that we have two different ways of generating the odd order magic square produced by magic(n). 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.
Our new algorithm generates the same magic square with a few array operations.
[I,J] = ndgrid(1:n);
A = mod(I+J+(n-3)/2,n);
B = mod(I+2*J-2,n);
M = n*A + B + 1;
It is not hard to prove that A and B have full rank and so it is likely that their sum has full rank as well.
There 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.
Does 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.
Doubly Even Order
For 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.
Let's look at the generation of magic(4). Start with the integers 1:16 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.
1* 2 3 4*
5 6* 7* 8
9 10* 11* 12
13* 14 15 16*
Now regard the marked elements as a one-dimensional array and flip that array end-for-end. This produces magic(4).
16* 2 3 13*
5 11* 10* 8
9 7* 6* 12
4* 14 15 1*
It 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.
Order Four
What 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 Numerical Recipes. Did Dürer Intentionally Show Only His Second-Best Magic Square? (I won't spoil Bill's punch line by revealing the answer to the question he poises in his title.)
In order to investigate Dürer'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, enumerate.m, that checks the rank of all magic squares of order 4. This is the first time I have ever written a MATLAB code with for loops nested seven deep. It generates
16*15*14*13*12*11*10
= 57657600
matrices. Of these, only 7040 pass the check in the inner block to qualify as legal magic squares.
Dividing 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
rank count
3 640
4 240
So, 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.
By the way, the entire computation by enumerate4 takes about 415 seconds on my Lenovo X201 2.67 GHz i7 laptop. The first nonsingular 4-by-4 magic square that enumerate4 finds is
1 2 16 15
13 14 4 3
12 7 9 6
8 11 5 10
Singly Even
The bar chart above shows that the rank for singly even magic squares, where n 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.
The first step of the singly even algorithm is to stack together four copies of the square of odd order $n/2$.
p = n/2;
A = magic(p);
M = [A A+2*p^2; A+3*p^2 A+p^2];
This M has equal column sums, but not equal row sums. And, since A has rank $n/2$, this M has rank $n/2+1$
The next step interchanges some blocks of elements to correct the row sums.
i = (1:p)';
k = (n-2)/4;
j = [1:k (n-k+2):n];
M([i; i+p],j) = M([i+p; i],j);
It turns out that this step does not change the rank. Now this M has equal row and column sums, but its diagonal sums are not quite right.
The final step swaps two pairs of elements to correct the diagonal.
评论
要发表评论,请点击 此处 登录到您的 MathWorks 帐户或创建一个新帐户。