Mercurial > hg > octave-lyh
annotate scripts/special-matrix/magic.m @ 9245:16f53d29049f
update copyright notices
author | John W. Eaton <jwe@octave.org> |
---|---|
date | Fri, 22 May 2009 10:46:00 -0400 |
parents | 853f96e8008f |
children | 1740012184f9 |
rev | line source |
---|---|
9245 | 1 ## Copyright (C) 1999, 2000, 2001, 2006, 2007, 2009 Paul Kienzle |
5827 | 2 ## |
3 ## This file is part of Octave. | |
4 ## | |
5 ## Octave is free software; you can redistribute it and/or modify it | |
6 ## under the terms of the GNU General Public License as published by | |
7016 | 7 ## the Free Software Foundation; either version 3 of the License, or (at |
8 ## your option) any later version. | |
5827 | 9 ## |
10 ## Octave is distributed in the hope that it will be useful, but | |
11 ## WITHOUT ANY WARRANTY; without even the implied warranty of | |
12 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU | |
13 ## General Public License for more details. | |
14 ## | |
15 ## You should have received a copy of the GNU General Public License | |
7016 | 16 ## along with Octave; see the file COPYING. If not, see |
17 ## <http://www.gnu.org/licenses/>. | |
5827 | 18 |
19 ## -*- texinfo -*- | |
20 ## @deftypefn {Function File} {} magic (@var{n}) | |
21 ## | |
9041
853f96e8008f
Cleanup documentation file matrix.texi
Rik <rdrider0-list@yahoo.com>
parents:
7017
diff
changeset
|
22 ## Create an @var{n}-by-@var{n} magic square. Note that @code{magic |
5827 | 23 ## (@var{2})} is undefined since there is no 2-by-2 magic square. |
24 ## | |
25 ## @end deftypefn | |
26 | |
27 function A = magic(n) | |
28 | |
29 if (nargin != 1) | |
30 print_usage (); | |
31 endif | |
32 | |
33 if (n != floor (n) || n < 0 || n == 2) | |
34 error ("magic: n must be an positive integer not equal to 2"); | |
35 endif | |
36 | |
37 if (n == 0) | |
38 | |
39 A = []; | |
40 | |
41 elseif (mod (n, 2) == 1) | |
42 | |
43 shift = floor ((0:n*n-1)/n); | |
44 c = mod ([1:n*n] - shift + (n-3)/2, n); | |
45 r = mod ([n*n:-1:1] + 2*shift, n); | |
46 A (c*n+r+1) = 1:n*n; | |
47 A = reshape (A, n, n); | |
48 | |
49 elseif (mod (n, 4) == 0) | |
50 | |
51 A = reshape (1:n*n, n, n)'; | |
52 I = [1:4:n, 4:4:n]; | |
53 J = fliplr (I); | |
54 A(I,I) = A(J,J); | |
55 I = [2:4:n, 3:4:n]; | |
56 J = fliplr (I); | |
57 A(I,I) = A(J,J); | |
58 | |
59 elseif (mod (n, 4) == 2) | |
60 | |
61 m = n/2; | |
62 A = magic (m); | |
63 A = [A, A+2*m*m; A+3*m*m, A+m*m]; | |
64 k = (m-1)/2; | |
65 if (k>1) | |
66 I = 1:m; | |
67 J = [2:k, n-k+2:n]; | |
68 A([I,I+m],J) = A([I+m,I],J); | |
69 endif | |
70 I = [1:k, k+2:m]; | |
71 A([I,I+m],1) = A([I+m,I],1); | |
72 I = k + 1; | |
73 A([I,I+m],I) = A([I+m,I],I); | |
74 | |
75 endif | |
76 | |
77 endfunction | |
78 | |
79 %!test | |
80 %! for i=3:30 | |
81 %! A=magic(i); | |
82 %! assert(norm(diff([sum(diag(A)),sum(diag(flipud(A))),sum(A),sum(A')])),0) | |
83 %! endfor | |
84 %!assert(isempty(magic(0))); | |
85 %!assert(magic(1),1); | |
86 %!error magic(2) |