Academic Articles & Supplements

IMfinal30Jan2021.nb

Peter D. Loly, Adam Rogers, Ian D. Cameron

Author

Peter D. Loly, Adam Rogers, Ian D. Cameron

Title

IMfinal30Jan2021.nb

Description

Supplemental notebook to "Frierson's 1907 Parameterization of Compound Magic Squares Extended to Orders 3^{l}, l=1,2,3,..., with Information Entropy"

Integer Matrices - The Science

Peter D. Loly

Singular Values (SVs) of order n square integer matrices are the basis of a spectral code for scientific study of Latin, magic and other integer squares.
This is critical because some squares may only have one or two ** non-zero eigenvalues (EVs), but as many SVs as their matrix rank.
TestDA checks the row and column line-sums for the doubly affine (DA) nature of these sums for both Latin squares and for magic squares, the latter with the same line-sums for the principal (d1) and counter (dexter) (d2) diagonals.
Pandiag tests whether all parallel broken diagonals (pandiagonals) have the same sum as their rows and columns.
The main code, Spectra, focusses on their SVs which lead to important classification measures using Shannon information entropy, and a related % measure for comparisons.
Finally code Compound combines smaller squares of orders m,n into larger ones of composite (product) order mn, whose properties can then be assessed with the earlier codes.
Also an order 4 magic cube is unfolded to an order 8 magic square, while a magic rectangle and a random (integer) square are also discussed.
The main Spectra code may also have some use for non-integer matrices, and even ones with complex elements, e.g. quantum angular momentum matrices.

Introduction to Integer Squares

There are many Wolfram-Mathematica links to magic squares (MSs) which will be referenced, including my 2007 publication in Complex Systems [1] of a presentation at the 2006 “New Kind of Science” [NKS] conference, which built on a 2003 conference talk [2] in the International History and Philosophy of Science Teaching (IHPST) series (Munich 1998, 2000, Winnipeg 2003). There I laid out the matrix eigenvalue (EV) and eigenfunction properties of magic squares, where the row-column linesum eigenvalue is both positive and largest, but others may be negative, complex or even vanish. This was followed by more detailed reports for the 880 set of order 4 and larger magic squares at IWMS Windsor 2007 in (LAA [3] Tables 5,6,7 (2009), which included matrix Singular Value properties (SVD - Singular Value Decomposition) motivated by MSs with just one non-zero EV, which was extended at Bedlewo (Poland) 2012 (as a video presentation) in DMPS 2013 [4], with the new tools now forming the basis of my main Spectra code in the present notebook, and with many tables for squares from orders 2 to 9.

The unpublished dataset from Bedlewo, in Mathematica’s bracket notation, is now made available as a Notebook [5].

I discuss the simpler Latin squares of Euler first:

Latin squares

Latinsquares[6]ofordernuseoneofeachsymbolineachrowandcolumn,herethesequenceofpositiveintegers,1,2,..n,withthesmallest:

In[]:=

latin2={{1,2},{2,1}};MatrixForm[latin2]

Out[]//MatrixForm=



1	2
2	1



whichhasrowandcolumnsumsof3,andwherethemaindiagonald1sumsto2,whilethecounter(ordexter)diagonald2sumsto4.GenerallytheRowColumn[RC]linesums(

)oftheseLatinsquaresaregivenby:

In[]:=

LatinRC[n_]:=N[n(n+1)/2]

for which LatinRC[2]=3.

There is also a single order 3 Latin square with an RC linesum of 6:

In[]:=

latin3={{1,2,3},{3,1,2},{2,3,1}};

Completed solutions of Sudoku puzzles are order 9 Latin squares partitioned into order 3 subsquares. Later I will show a simple order 4 Latin square with order 2 subsquares which I call mini-Sudoku, sud4a.

Magic squares

Magicsquares[7]usingthesequenceofpositiveintegers,1,2,..,

,beginwiththesoleorderthreeancientmagicsquare,theLuoshu,documented~369-286B.C.E.[8](Swetz)withitsverticalinvertmagic3:

In[]:=

magic3={{8,1,6},{3,5,7},{4,9,2}};MatrixForm[magic3]

Out[]//MatrixForm=

8	1	6
3	5	7
4	9	2

where the main diagonal d1=8+5+2=15 and the counter diagonal d2=6+5+4=15.

The Row - Column - Diagonal [RCD] linesums (

) of magic squares are given by:

In[]:=

MagicRCD[n_]:=N[n(

+1)/2]

which for magic3 gives MagicRCD[3]=15.

A search “wolfram magic square” found MagicSquare.nb by Singh [9] (20 February 2019) for a magic square of dimension n. I also found a magic[n] macro for a single magic square of any odd, even and double even order in Mathematica Stack Exchange by Nasser[10] (2015), together with an improvement by "chyanog", and gives magic3 above.

Background

An early connection to magic square matrices was made in a short note in 1877 by Arthur Cayley [11], one of the founders of matrix algebra. At order n=4 Frénicle (1605-1675, published in 1693!) found that there are 880 distinct magic squares with RCD’s of 34, which Dudeney [12] classified into 12 groups c. 1900 - see Andrews classic text [13]. Then for a complete list using the elements 1,2 ,...16 see Benson & Jacoby [14], Heinz [15] and Heinz & Hendricks ]16], amongst them the famous magic square of Albrecht Dürer, in which the middle elements of the last row celebrate 1514:

In[]:=

durer={{16,3,2,13},{5,10,11,8},{9,6,7,12},{4,15,14,1}};

For order 5 there are some 275 million distinct magic squares - see Walter Trump [17] for estimates of populations of higher orders.

Doubly-Affine [DA] square matrices - TestDA

These have the same row and column sums. Magic squares (and diagonal Latin squares) in addition to being DA, have the same linesum for both main diagonals. Often these properties are known beforehand. Clearly small squares can be quickly assessed for row, column and diagonal line-sums by visual inspection, but larger squares and large numbers of a given type need systematic code for checking the characterising properties.
This code assumes a square matrix and calculates its row, column and main diagonal line-sums using standard matrix and vector operations:

In[]:=

TestDA[mat_] := (fi = mat; dim = Dimensions[fi]; n = dim[[1]];m = dim[[2]]; en=ConstantArray[1,{n,1}];rowsums = fi.en;em=ConstantArray[1,{1,n}];colsums = em.fi;d1=0.;For[i = 1, i < n + 1, i++, d1 = d1 + fi[[i, i]]]; d2=0.;For[j = 1, j < n + 1, j++, d2 = d2 + fi[[n + 1 - j, j]]]; s=MagicRCD[n];l=LatinRC[n];Print[fi//MatrixForm,", d1=",d1,", d2=",d2", rows=",rowsums,", cols=",colsums]);

Since there is just a single Latin square of orders 2 and 3 (aside from rotations and reflections), these are tested first:

In[]:=

TestDA[latin2]



1	2
2	1

, d1=2., d2=4., rows={{3},{3}}, cols={{3,3}}

Not printed in the above for brevity of output are:

=MagicRCD[n];

=LatinRC[n], which can be accessed with:

In[]:=

Print["MagicRCD=",s,", LatinRC=",l]

MagicRCD=5., LatinRC=3.

Next the 3rd order Latin square:

In[]:=

TestDA[latin3]

1	2	3
3	1	2
2	3	1

, d1=3., d2=6., rows={{6},{6},{6}}, cols={{6,6,6}}

Note that the diagonal sums are different from the RC sums in both cases. Rotation and/or reflection of latin2,3 just swaps d1 and d2.

In[]:=

TestDA[magic3]

8	1	6
3	5	7
4	9	2

, d1=15., d2=15., rows={{15},{15},{15}}, cols={{15,15,15}}

now with both diagonals equal to

=15 as per MagicRCD given earlier.

Integer Sequence Doubly-Affine matrices [ISDA]

Magic and Latin squares are examples of Integer Sequence Doubly-Affine (ISDA’s) matrices - see (RCL)[18], which used 0,1,2,..,

(or

. Here they will be taken to run 1,2,..,

(or

Pandiagonality - code Pandiag

Also called panmagic. To check this duplicate a square matrix vertically (or horizontally) to a rectangular one, e.g. :

In[]:=

vertpair=Join[magic3,magic3];MatrixForm[vertpair]

Out[]//MatrixForm=

8	1	6
3	5	7
4	9	2
8	1	6
3	5	7
4	9	2

Here the pandiagonals parallel to d1=8+5+2=15 are d1r=1+7+4=12 and d1rr=6+3+9=18, and those parallel to d2=6+5+4=15 are d2l=7+9+8=24 and d2ll=2+1+3=6, which for both cases average 15.
My code to check if square matrices have the same pandagonal linesums:

In[]:=

Pandiag[mat_]:=(fi=mat;dim=Dimensions[fi];xfi=Join[fi,fi];n=dim[[2]];check=1;d1=0.;For[i=1,i<n+1,i++,d1=d1+xfi[[i,i]]];pd=0.;For[i=1,i<n+1,i++,pd=pd+xfi[[i+n-1,i]]];If[pd≠d1,check=0];dp=0.;For[j=1,j<n+1,j++,dp=dp+xfi[[n+1-j,j]]];If[dp≠d1,check=0];Print["check=",check]);

where check = 1 will indicate pandiagonal, and check = 0 not pandiagonal.

N.B. “vertpair” effects the periodic continuation of the matrix in the vertical direction - an idea used in crystal physics theory to eliminate end effects. Alternatively one may continue in the other direction:

In[]:=

horizpair={{8,1,6,8,1,6},{3,5,7,3,5,7},{4,9,2,4,9,2}}//MatrixForm

Out[]//MatrixForm=

8	1	6	8	1	6
3	5	7	3	5	7
4	9	2	4	9	2

Some examples

For magic4 [9] which is Frénicle F176 [14] :

In[]:=

magic4={{16,2,3,13},{5,11,10,8},{9,7,6,12},{4,14,15,1}};

In[]:=

Pandiag[magic4]

check=0

and so not pandiagonal, whereas Weisstein’s example [19] (Frénicle F102 Dudeney Group I [13]) is:

In[]:=

F102={{1,8,10,15},{12,13,3,6},{7,2,16,9},{14,11,5,4}};

In[]:=

Pandiag[F102]

check=1

Dürer’s magic square given earlier is not pandiagonal, nor are magic3,4 pandiagonal, although 48 of the 880 order 4’s are pandiagonal [3,19].

Torus

Note that the rectangles duplicating a matrix may be considered equivalent to wrapping around either radius of a torus (bagel) [20]. This is tantamount to periodic boundary conditions which I first encountered in connection in the context of calculations for crystals, which would for n=4 wrap around the smaller radius, and that repeated around the larger radius:

In[]:=

ParametricPlot3D[{(2+Cos[v])Cos[u],(2+Cos[v])Sin[u],Sin[v]},{u,0,2Pi},{v,0,2Pi},Mesh3,BoundaryStyleBlack,PlotStyleFaceForm[Red,Yellow]]

Out[]=

Here Mesh=3 gives the wrapping needed for an order 4 square - use Mesh=2 for n=3, and Mesh=7 for chessboard squares of n=8, etc.

Associative magic squares

Another important type - the associative magic squares in which all pairs antipodal about the centre all have the sum 1+

[21], which are quite easy to check visually, and for orders 3,4,5,6
number 1, 48, 48544, 0 respectively. Ultramagic squares are both pandiagonal and associative.

Most-perfect pandiagonal (MPPD) magic squares

Introduced by McClintock [22] in 1897, and the definitive study by Ollerenshaw and Brée [23] in 1998. They used sequences 0,1,2,...,

-1, so I add 1 to their elements to conform to the more common 1,2,..,

used in the majority of references herein.
MPPDs have the additional property that:
(i) all 2-by-2 sub-squares have the same sum, including those that run over the edges when tiled or when wrapped over a torus, and
(ii) each element is complementary to the one distant from it two steps along any pandiagonal, e.g. from Harry White [24]:

In[]:=

mppd4={{1,8,10,15},{12,13,3,6},{7,2,16,9},{14,11,5,4}};MatrixForm[mppd4]

Out[]//MatrixForm=

1	8	10	15
12	13	3	6
7	2	16	9
14	11	5	4

and Pandiag[mppd4] gives check=1.

From [25]: “The enumeration of all ‘most-perfect’ squares of order n = 2r, (r > 1) that are pandiagonal magic squares with additional special properties and are defined in Chapter 1, was conjectured by me (Ollerenshaw) in 1987 and stated in the 25th Anniversary IMA Bulletin of March 1989. Formal publication of the proof was delayed because extensions of this result for other multiples of 4 were in sight. The full enumeration for n =

(p any prime > 2, r > 1, s ≥ 0) that forms the first part of this book was arrived at in 1989 and its existence mentioned in the IMA Bulletin article referred to above.” Which includes orders 8, 12, 20,... etc.N.B. so far PDL has only found these important MSs on Wolfram Mathematica sites in my NKS2006 paper [1].

Franklin’s bent diagonal magic squares on a Chessboard (n=8)

Students Schindel and Rempel with Loly c. 2004/5 made an (exact) backtracking count of Ben Franklin’s bent diagonal magic squares on a chessboard (1750), exactly three times that of the order 8 MPPDs (368,640), which we published in PRSA [26] in 2006 on the 300th anniversary of Franklin’s birth, and was featured in my later NKS2006 talk [1] published in 2007. See later for a more recent backtracking program by Robert Cowen in TMJ[27]).

Franklin’s bent diagonal magic squares on a Chessboard (n=8)

Other integer squares as matrices

While our main interest concerns integer magic squares, the simpler integer Latin squares, as well as integer square addition and multiplication tables, our main matrix code, Spectra, has some use for other matrices, including random integer squares and rectangular matrices.

Some basic matrices

Constant square matrices can be obtained from:

In[]:=

csq[n_,c_]:=ConstantArray[c,{n,n}]

In[]:=

csq[2,1]//MatrixForm

Out[]//MatrixForm=



1	1
1	1



Then the identity matrix by:

In[]:=

Id[n_]:=IdentityMatrix[n]

In[]:=

Id[2]//MatrixForm

Out[]//MatrixForm=



1	0
0	1



and its rotation, left-right or up-down flip obtained by using:

In[]:=

Jsq[n_]:=Reverse[Id[n]]

In[]:=

Jsq[2]//MatrixForm

Out[]//MatrixForm=



0	1
1	0



To see how the flip works take a general 2-by-2 matrix:

In[]:=

abcd={{a,b},{c,d}};

first a left-right flip:

In[]:=

Dot[abcd,Jsq[2]]//MatrixForm

Out[]//MatrixForm=

b	a
d	c

and then a top-bottom flip:

In[]:=

Dot[Jsq[2],abcd]//MatrixForm

Out[]//MatrixForm=

c	d
a	b

Square Addition Tables

Addition tables where the successive rows are continued increments of the previous row, e.g. the smallest:

In[]:=

at2={{1,2},{3,4}};MatrixForm[at2]

Out[]//MatrixForm=



1	2
3	4



at2 is pandiagonal with d1=d2=5, but with no rows or columns with that sum.

Next we show the pandiagonality of the third order square addition table, at3, which is the numerical keypad on telephone dial pads, as well as its vertical reflection on computer keyboards:

In[]:=

at3={{1,2,3},{4,5,6},{7,8,9}};MatrixForm[at3]

Out[]//MatrixForm=

1	2	3
4	5	6
7	8	9

which checks to be pandiagonal, as are all larger addition tables.

Beyond these addition tables our attention was drawn to the next types in 2019 for which we next study their matrix properties:

More integer squares - Sum and Product Grids

Recently Hartnett [28] discussed Sum Grids and Product Grids, the smallest two Sum Grids show their trend:

In[]:=

sum2={{3,4},{2,3}};MatrixForm[sum2]

Out[]//MatrixForm=



3	4
2	3



In[]:=

sum3={{4,5,6},{3,4,5},{2,3,4}};MatrixForm[sum3]

Out[]//MatrixForm=

4	5	6
3	4	5
2	3	4

both also pandiagonal.
By contrast the first Product Grid is not pandiagonal (obvious by inspection!):

In[]:=

prod2={{2,4},{1,2}};MatrixForm[prod2]

Out[]//MatrixForm=



2	4
1	2



nor is the next:

In[]:=

prod3={{3,6,9},{2,4,6},{1,2,3}};MatrixForm[prod3]

Out[]//MatrixForm=

3	6	9
2	4	6
1	2	3

Matrix spectra - Eigenvalues (EVs) and Singular Values (SVs)

Our studies the EV spectra of magic squares deepened at the 2007 IWMS-16 conference [3] (LAA) by including Singular Value (SV) spectra whose squares are the eigenvalues of the Gramian product of the matrix and its transpose, and which provide a set of decreasing positive SVs whose number equals the matrix rank, r. The SVs are often more useful than the eigenvalues, especially when some or most EVs vanish - more on this later. A critical leap followed from reading in 2010 the forthcoming abstract by Newton and DeSalvo (NDS)[29] who studied the Shannon [30] entropy of Sudoku matrices, in Proc. R. Soc. A. (PRSA), whose abstract I read as a result of our earlier paper [26] in the same journal - Sudoku’s being order 9 Latin squares caught my attention and it was easy to use their measures based on Claude Shannon’s information entropy for other Latin and magic squares, as well as now for other integer matrices. A useful % measure called Compression, C, bounded between 0% and 100%, facilitates comparisons between different squares. See our data [5].After reading background in Miller’s history [31] I saw how to extend NDS’s measures by adapting formulae of Albert Girard in 1629 for the sums of powers of the coefficients of the eigenvalue characteristic polynomial, which are integer for integer squares, to the coefficients of the singular values characteristic polynomial to give sums of their even powers. This was presented in a video sent to 2012 conference at Bedlewo, Poland and published in the proceedings DMPS [4], but that publisher omitted my advertised data appendix of all the matrices used therein (in Mathematica’s bracket notation), which formed the basis of extensive tables. That data has now been put in Notebook style [5] to provide readers with the input matrices used in [4].

Code Spectra

This also uses equations first derived in DMPS [4] for the coefficients b1,b2,b3 of the characteristic equation of the Gramian matrix, which is the product of a matrix and its transpose. For an integer matrix, the b1,b2,b3 coefficients of Gramian are also integers and lead to expressions for the sums of the 4th power of the SVs, L, which I found also to be integer for integer squares - however L can grow quite large!
L’s reduction for DA (Doubly Affine) squares to the shorter integer R which omits the large 4th power of the first (RC linesum) SV, and also Q, the sum of the 6th powers of the SVs less its first, which was used once in DMPS [4] to resolve a degeneracy of L (and R) in the order 4 set of 880 distinct magic squares. See my code below for formulae for L,R,Q.
Note that if the series begins with 0,1,.. instead of 1,2,.. then L changes but R,Q do not, but of course the Compression, C, and the entropy H would change.
The full set of matrix elements for the 63 SV clans with their Frénicle index numbers and Dudeney Groups includes clans 25 and 26 which have the same R value but different ranks, which can be distinguished by Q described above.
Included in this spectra code are two related measures inspired by the famous entropy equation of Boltzmann on his grave in Vienna [32],

log(W)

, here from the entropy which NDS labelled “ht”, and labelled “H” in my code.
Also included is a simple “Spread” measure of the matrix elements which is normalized by the linesum EV, ev[[1]], which only fails in a rare event if that EV vanishes!

In[]:=

Spectra[mat_]:=(f=mat;dim=Dimensions[f];n=dim[[1]];cpev=CharacteristicPolynomial[f,x];r=MatrixRank[f];ev=Eigenvalues[f];gram=f.Transpose[f];svsqd=Eigenvalues[gram];cpSigmaSq=CharacteristicPolynomial[gram,x];b1=-(-1)^n*Coefficient[cpSigmaSq,x,n-1];b2=(-1)^n*Coefficient[cpSigmaSq,x,n-2];b3=-(-1)^n*Coefficient[cpSigmaSq,x,n-3];L=b1^2-2*b2;sigma=SingularValueList[f];SingularValueDecomposition[f];R=L-sigma[[1]]^4;Q=b1^3-3*b1*b2+3*b3-sigma[[1]]^6;sigmaTot=0.;For[p=1,p<r+1,p++,sigmaTot=sigmaTot+sigma[[p]]];sigmaN=sigma/sigmaTot;H=0.;For[j=1,j<r+1,j++,H=H-sigmaN[[j]]*Log[sigmaN[[j]]]];comp=(1-H/Log[n])*100;spread=n*(Max[f]-Min[f])/ev[[1]]//N;Print["n=",n,", rank=",r,", C=",comp,"%, H=",H];);

The printout of spectra focuses on order n, rank r, Compression C %, and entropy H, first for the first and oldest magic square :

In[]:=

Spectra[magic3]

n=3, rank=3, C=14.7017%, H=0.937098

where a compact default printout gives 6 digit for both C and H. This can be supplemented by statements needed to print other ingredients as desired, e.g. Print[ev or Nev, sv or sigmas, ...], e.g.

In[]:=

Print["EV=",N[ev],", charpolyEV:",cpev,", charpolySV:",cpSigmaSq,", SVsqd=",N[svsqd],", L=",L,", R=",R,", Q=",Q,", Spread=",spread]

EV={15.,-4.89898,4.89898}, charpolyEV:-360+24x+15

, charpolySV:129600-14076x+285

, SVsqd={225.,48.,12.}, L=53073, R=2448, Q=112320, Spread=1.6

where summing the squares of the squared SVs explicitly:

(225)

(48)

(12)

=53073, or the shorter: R=L-

(225)

=2448, and Q=

(48)

(12)

=112320.Here the Spread is calculated from 3(9-1=8)/15=8/5=1.6 (exactly).For magic4, which is a flip about the counter diagonal of F176 [15], has same spectra as Frénicle 21, Dudeney Group 2, and has SV clan 9 [4,6]:

Spectra[magic4]

n=4, rank=3, C=37.2284%, H=0.8702

where I note the minimum rank 3 of magic squares which Drury [33] obtained after I gave a talk at McGill University.

Print[ev,svsqd]

34,-4

,0{1156,320,20,0}

For the 880 order 4 magic squares there are 63 SV clans (of which DMPS p.127, Table 1 showed only 13 of these clans).N.B. Different values for C%,H will be found if the squares do not use elements 1,2,3,...Next a few select examples will now be presented to show the scope and power of these codes.

A magic square with a single non-zero EV - “1EV”

One of the most interesting order 4 magic squares is F360 (Dudeney [13] Group 3 and SV clan 9 of the 63 clans:

In[]:=

F360={{2,11,7,14},{13,8,12,1},{16,5,9,4},{3,10,6,15}};MatrixForm[F360]

Out[]//MatrixForm=

2	11	7	14
13	8	12	1
16	5	9	4
3	10	6	15

which is clearly associative in that all antipodal elements about the centre sum to 17. Then:

In[]:=

Spectra[F360]

n=4, rank=3, C=37.2284%, H=0.8702

and for more details:

In[]:=

Print[ev,svsqd]

{34,0,0,0}{1156,320,20,0}

showing just one EV, but the minimum rank 3 of magic squares [33] as shown by the SVs, where 1156 is 34-squared. Note the positive decreasing SVs whose number of non-zero values gives the matrix rank, here 3.

Next the first two integer squares of all ones:

In[]:=

e2=csq[2,1]={{1,1},{1,1}};e3=csq[3,1];e4=csq[4,1];

and all orders of these have full 100% Compression, e.g.:

In[]:=

Spectra[e2]

n=2, rank=1, C=100.%, H=0.

This 100% Compression represents complete order, i.e. all elements identical, therefore the lowest entropy H=0, and is true for all larger squares of all 1’s. Then all other non-constant matrices have lower order and higher entropy.

In[]:=

Print[ev,svsqd]

{2,0}{4,0}

Then all addition and sum tables are rank 2 with two non-zero eigenvalues and singular values - “2EV” and “2SV”, e.g.:

In[]:=

Spectra[at2]

n=2, rank=2, C=66.1685%, H=0.234502

In[]:=

Print[ev,svsqd]



(5+

(5-

)15+

221

,15-

221



Next my mini-Sudoku Latin square with just one eigenvalue [1EV] but three SVs (rank 3) [4]:

In[]:=

sud4a={{1,2,3,4},{3,4,1,2},{4,3,2,1},{2,1,4,3}};

In[]:=

TestDA[sud4a]

1	2	3	4
3	4	1	2
4	3	2	1
2	1	4	3

, d1=10., d2=10., rows={{10},{10},{10},{10}}, cols={{10,10,10,10}}

which is seen to be RCD Latin (and “semi-pandiagonal” in alternate ones), and which I termed “mini-Sudoku” since 1,2,3,4 occur in each row, column and each tiled 2-by-2 subsquare [4], hence sud4a.

In[]:=

Pandiag[sud4a]

check=0

but can be considered as semi-pandiagonal since 2nd, 4th, ... (alternate) pandiagonals have the RC sum!

In[]:=

Spectra[sud4a]

n=4, rank=3, C=35.0603%, H=0.900256

In[]:=

Print[ev,svsqd]

{10,0,0,0}{100,16,4,0}

showing just one EV, but three non-zero SVs.

Sum and Product squares

Hartnett [28] begins with sum2 and prod2, shown earlier and now having:

In[]:=

Spectra[sum2]

n=2, rank=2, C=82.7872%, H=0.11931

In[]:=

Print[ev,svsqd]

3+2

,3-2

19+6

19+6



with quite different spectra for:

In[]:=

Spectra[prod2]

n=2, rank=1, C=100.%, H=0.

In[]:=

Print[ev,svsqd]

{4,0}{25,0}

and with rank=1 and 100% Compression continuing for higher orders.

Then:

In[]:=

Spectra[sum3]

n=3, rank=2, C=85.5652%, H=0.158583

In[]:=

Print[ev,svsqd]

6+

,6-

,06(13+2

13+2

,0

showing an increased Compression for the “sums” as n increases - NEV!

Compounding pairs of Latin or magic squares to larger product orders - GAPDA? RCL

A new world record for the largest magic square of order 12,544 was created in 2002 with Chan [34].
I used the present Compound code in my calculations for RCL [18]. It is intended to be used for pairs of Latin squares with k=1, or pairs of magic squares with k=2, but never a mixture!

In[]:=

Compound[ik_,mat1_,mat2_]:=(k=ik;sq1=mat1;sq2=mat2;dimm=Dimensions[sq1];dimn=Dimensions[sq2];m=dimm[[1]];n=dimn[[1]];Print["n=",n,", m=",m," mat1=",mat1," mat2=",mat2];alpha=KroneckerProduct[sq1,csq[n,1]];beta=KroneckerProduct[csq[m,1],sq2];cA=(n^k)*alpha+beta-(n^k)*csq[n*m,1];cD=alpha+(m^k)*beta-(m^k)*csq[n*m,1];);

Compound calculates a pair of order mn product squares called cA,cD (which can then be accessed by those names) - the first called “Associated” as the sub-squares are intact and only increased in overall magnitude:

In[]:=

Compound[2,magic3,magic3]

n=3, m=3 mat1={{8,1,6},{3,5,7},{4,9,2}} mat2={{8,1,6},{3,5,7},{4,9,2}}

cA here is a very ancient magic square - see Swetz [9]:

In[]:=

MatrixForm[cA]

Out[]//MatrixForm=

71	64	69	8	1	6	53	46	51
66	68	70	3	5	7	48	50	52
67	72	65	4	9	2	49	54	47
26	19	24	44	37	42	62	55	60
21	23	25	39	41	43	57	59	61
22	27	20	40	45	38	58	63	56
35	28	33	80	73	78	17	10	15
30	32	34	75	77	79	12	14	16
31	36	29	76	81	74	13	18	11

where the incremented sub-squares of mat1 stand out clearly. Then its partner cD has the elements of mat1 Distributed over different sub-squares of mat2:

In[]:=

MatrixForm[cD]

Out[]//MatrixForm=

71	8	53	64	1	46	69	6	51
26	44	62	19	37	55	24	42	60
35	80	17	28	73	10	33	78	15
66	3	48	68	5	50	70	7	52
21	39	57	23	41	59	25	43	61
30	75	12	32	77	14	34	79	16
67	4	49	72	9	54	65	2	47
22	40	58	27	45	63	20	38	56
31	76	13	36	81	18	29	74	11

cD has the same spectra as cA:

In[]:=

Spectra[cA]

n=9, rank=5, C=48.572%, H=1.12999

In[]:=

Print[ev,svsqd]

369,-54

,54

,-6

,0,0,0,0{136161,34992,8748,432,108,0,0,0,0}

These squares were both given in 1907 by Frierson in Andrews [14], along with two other pairs, which recently Loly and Cameron [36] showed to have smaller Compressions of 40.0241% and 39.8296% - i.e. somewhat less ordered.
A general result was found by Rogers, Cameron and Loly (RCL)[35] that the ranks of compound squares are the sum of their component ranks less 1, here 3+3-1=5. This also applies to repeated compounding [35].

Compounding mppd4

Now return to mppd4 and first check its spectra before compounding:

In[]:=

Spectra[mppd4]

n=4, rank=3, C=30.8812%, H=0.95819

In[]:=

Compound[2,mppd4,mppd4]

n=4, m=4 mat1={{1,8,10,15},{12,13,3,6},{7,2,16,9},{14,11,5,4}} mat2={{1,8,10,15},{12,13,3,6},{7,2,16,9},{14,11,5,4}}

and Pandiag[cA] gives check = 1, so cA is pandiagonal.

In[]:=

Spectra[cA]

n=16, rank=5, C=60.8479%, H=1.08553

So the Compression nearly doubles from 30.8812% to 60.8479%

Frierson used an algebraic approach to compounding order 3 magic squares

While focussed on integer squares it is worth noting that compounding can be used algebraically. See our recent work on Frierson [36] for his algebra which begins with order 3:

In[]:=

fr3[v,y]:={{2v+y,0,v+2y},{2y,v+y,2v},{v,2y+2v,y}};

In[]:=

MatrixForm[fr3[v,y]]

Out[]//MatrixForm=

2v+y	0	v+2y
2y	v+y	2v
v	2v+2y	y

In[]:=

2v+y	0	v+2y
2y	v+y	2v
v	2v+2y	y

Out[]=

{{2v+y,0,v+2y},{2y,v+y,2v},{v,2v+2y,y}}

In[]:=

ff3[ss,t]:={{2ss+t,0,ss+2t},{2t,ss+t,2ss},{ss,2t+2ss,t}};

In[]:=

MatrixForm[ff3[ss,t]]

Out[]//MatrixForm=

2ss+t	0	ss+2t
2t	ss+t	2ss
ss	2ss+2t	t

then we obtain Frierson’s algebra for order 9 using our compound code:

In[]:=

Compound[1,ff3[ss,t],fr3[v,y]]

n=3, m=3 mat1={{2ss+t,0,ss+2t},{2t,ss+t,2ss},{ss,2ss+2t,t}} mat2={{2v+y,0,v+2y},{2y,v+y,2v},{v,2v+2y,y}}

Our recent arXiv [35] has more detail, including extensions to higher powers of 3 which indicate a fractal limit.

Latin squares may also be compounded

In[]:=

Compound[1,latin2,latin2]

n=2, m=2 mat1={{1,2},{2,1}} mat2={{1,2},{2,1}}

which is easy to appreciate in its basic case:

In[]:=

MatrixForm[cA]

Out[]//MatrixForm=

1	2	3	4
2	1	4	3
3	4	1	2
4	3	2	1

In[]:=

Spectra[cA]

n=4, rank=3, C=35.0603%, H=0.900256

and while the same Compression as our earlier sud4a, it is not a mini-Sudoku, as it does not have one of each element in each 2-by-2 subsquare, although it does have the same spectra as lat4sud in DMPS [4].

Multimagic squares

The smallest bimagic square has order 8 and the smallest trimagic order 12 - see Eggermont [36], Derksen, Eggermont and van den Essen [37], Boyer [38], Weisstein [40], and van den Essen [41].
However we can explore the Latin square parallels with an order 5 Knut-Vik Latin square - see Nordgren [41], which remains pandiagonal when all its elements are raised to any power because there is always one of each in every row, column and pandiagonal:

In[]:=

KV5={{1,2,3,4,5},{4,5,1,2,3},{2,3,4,5,1},{5,1,2,3,4},{3,4,5,1,2}};

In[]:=

Pandiag[KV5]

check=1

testDA shows that this has RCDs of 15, and checks out also to be pandiagonal.
Then our spectra code gives:

In[]:=

Spectra[KV5]

n=5, rank=5, C=16.6093%, H=1.34212

I gave a conference talk on this at a CMS conference [42] in 2014, and leave the reader to check that compounding to order 25 preserves pandiagonality.
Now square all its elements:

In[]:=

KV5sq={{1,4,9,16,25},{16,25,1,4,9},{4,9,16,25,1},{25,1,4,9,16},{9,16,25,1,4}};

which satisfies our criteria for Latin squares to have one symbol in each row and column, and in fact has all RCDs of 55. Then this is still pandiagonal, with rank 5, but less than half the Compression at C=7.41065..%, H=1.49017.., spread=2.18182..

Then to show that spectra code can be used for non-integer squares take the square roots of the elements of kv5, i.e.:

In[]:=

r2=N[Sqrt[2]];r3=N[Sqrt[3]];r5=N[Sqrt[5]];

In[]:=

KV5sqrt={{1,r2,r3,2,r5},{2,r5,1,r2,r3},{r2,r3,2,r5,1},{r5,1,r2,r3,2},{r3,2,r5,1,r2}};

and find that the rank is still full at 5, but with almost double the Compression of KV5 at C=31.81834%, H=1.09734.

Backtracking

Earlier I noted a backtracking calculation for Franklin squares [26] and now note another in this journal for solving order 9 Sudoku puzzles by Cowen in The Mathematica Journal [27].

In[]:=

Cowen9S={{8,1,2,7,5,3,6,4,9},{9,4,3,6,8,2,1,7,5},{6,7,5,4,9,1,2,8,3},{1,5,4,2,3,7,8,9,6},{3,6,9,8,4,5,7,2,1},{2,8,7,1,6,9,5,3,4},{5,2,1,9,7,4,3,6,8},{4,3,8,5,2,6,9,1,7},{7,9,6,3,1,8,4,5,2}};

Once solved these can be analyzed with the codes presented above. I find for this one of Cowen' s full rank 9 Sudoku Latin square the expected linesum of

= 45.

In[]:=

Spectra[Cowen9S]

n=9, rank=9, C=23.4398%, H=1.6822

which C% is in the middle range of DMPS [4] Table 8.

Unfolding magic cubes to magic squares

Consider the order 4 magic cube from page 78 in Andrews [14]. This can be unfolded to an order 8 magic square:

In[]:=

Cox4={{1,63,62,4,48,18,19,45},{60,6,7,57,21,43,42,24},{56,10,11,53,25,39,38,28},{13,51,50,16,36,30,31,33},{32,34,35,29,49,15,14,52},{37,27,26,40,12,54,55,9},{41,23,22,44,8,58,59,5},{20,46,47,17,61,3,2,64}};

for which the RCDs are 260, n=8, rank=3, C=53.9919%, H=0.956711.., and appears to be associative, but is not pandiagonal.

This has a Compression falling about midway in the n=8 magic square Table 3 of DMPS[4].

Magic rectangles

I note that singular values of large rectangular integer matrices are also the basis of “Web” search engines - see also Langville [43]. Here I take a magic rectangle and instead of EVs form its two Gramians, one order 2 and the other order 4, both square matrices:

In[]:=

Rect={{8,2,3,5},{1,7,6,4}};sq22=Rect.Transpose[Rect];sq44=Transpose[Rect].Rect;

both have RCs of 162, sq22 has n=2, rank=2, C=26.6462%, H=0.956711, while sq44 has n=4, rank=2, C=63.3231%, H=0.50845.

Random Matrices

Following his publication of “New Kind of Science” (NKS) [44], Stephen Wolfram at the NKS conference in 2006 asked the present author about random squares - so here is a brief look at random integer squares:

In[]:=

Ran4=RandomInteger[16,{4,4}]

Out[]=

{{2,8,11,5},{3,7,9,12},{8,6,10,13},{0,16,9,4}}

In[]:=

TestDA[Ran4]

2	8	11	5
3	7	9	12
8	6	10	13
0	16	9	4

, d1=23., d2=20., rows={{26},{31},{37},{29}}, cols={{13,37,39,34}}

Ran4’s Compression is similar to Frénicle F46 in DMPS [4] Table1, but without any magic rows or columns.
For comparison with earlier magic squares add one to each element:

In[]:=

Spectra[Ran4+e4]

n=4, rank=4, C=33.2854%, H=0.924861

In[]:=

Print[N[ev],N[svsqd],", L=",L,", R=",R]

{35.5769,-9.53099,3.94037,-2.98627}{1323.58,137.574,13.2392,6.60385}, L=1771017, R=1771017-

1.32

…

where the lead EV is not an integer.
Caution - each time RandomInteger is activated the random integers change, and so also the spectra!

Conclusion

In reorganizing my codes for this paper I have created useful tools for exploring new issues, e.g. the recent Quanta Magazine article by Hartnett [28], as well as preparing to finish several projects presented at earlier conferences. While focussed on integer matrices, my codes can be used for others where some of my measures will no longer be integer, as indicated above.
A challenge remaining is to calculate the number of distinct SV clans for the 275 million order 5 magic squares! That would involve a backtracking calculation along the lines of Cowen’s [27] combined with the present tools.
In closing I note another integer context:

OEIS - Online Encyclopedia of Integer Sequences [45]

In an earlier study [46] I had obtained a formula for the physical moment of inertia of magic squares of any order n by replacing the matrix elements by corresponding amount of unit masses so that the mass of each row and column is the same, and then finding an exact formula for all orders, before another with Adam Rogers for magic cubes [47]. Jonathan Vos Pos filed those integer sequences with Neal Sloane’s OEIS - Online Encyclopedia of Integer Sequences [45] on Dec 23, 2006. Later Rogers and I found formulae for the electric multipole moments of the corresponding electrical problems [48].

Powers of 1EV squares

Earlier I showed that F360 and sud4a had just 1EV, in a recent arXiv paper with Cameron and Rogers [49] I found that powers of these quickly became constant matrices.

Another pandiagonal 4*4 square

An interest in the Myers-Briggs Type Table (MBTI) [50] interacted with periodic boundary conditions used to simplify crystal physics computations:

In[]:=

mbti={{8,6,2,4},{7,5,1,3},{15,13,9,11},{16,14,10,12}};

In[]:=

TestDA[mbti]

8	6	2	4
7	5	1	3
15	13	9	11
16	14	10	12

, d1=34., d2=34., rows={{20},{16},{48},{52}}, cols={{46,38,22,30}}

so not magic but is pandiagonal:

In[]:=

Pandiag[mbti]

check=1

In[]:=

Spectra[mbti]

n=4, rank=2, C=78.3964%, H=0.29949

Pauli matrices of quantum spin 1/2

I note that apart from a coefficient proportional to Planck’s constant [51], the Pauli spin 1/2 matrices of quantum mechanics, specifically its x-component:

In[]:=

sigX=latin2-e2;

In[]:=

MatrixForm[sigX]

Out[]//MatrixForm=



0	1
1	0



In[]:=

Spectra[sigX]

n=2, rank=2, C=0.%, H=0.693147

In[]:=

Print["EV=",N[ev],", charpolyEV:",cpev,", charpolySV:",cpSigmaSq,", SVsqd=",N[svsqd],", L=",L,", R=",R,", Q=",Q,", Spread=",spread]

EV={-1.,1.}, charpolyEV:-1+

, charpolySV:1-2x+

, SVsqd={1.,1.}, L=2, R=1, Q=1, Spread=-2.

For the record the other Cartesian components are

={{0,-i},{i,0}}

, where

is the square root of “-1”, and

={{1,0},{0,-1}}

which signals the up and down states. Also

iσ

=2{{0,1},{0,0}}

and

iσ

=2{{0,0},{1,0}}

I have studied higher spin angular momentum matrices which have larger dimensions but not pursued here: 2016 WCLAM 2016: “Quantum Angular Momentum Matrices” - see [51].

In[]:=

sigZ={{1,0},{0,-1}};Spectra[sigZ]

n=2, rank=2, C=0.%, H=0.693147

In[]:=

Print["EV=",N[ev],", charpolyEV:",cpev,", charpolySV:",cpSigmaSq,", SVsqd=",N[svsqd],", L=",L,", R=",R,", Q=",Q,", Spread=",spread]

EV={-1.,1.}, charpolyEV:-1+

, charpolySV:1-2x+

, SVsqd={1.,1.}, L=2, R=1, Q=1, Spread=-4.

In[]:=

sigPlus=2{{0,1},{0,0}};Spectra[sigPlus]

Power

:Infinite expression

encountered.

n=2, rank=1, C=100.%, H=0.

where the error message is caused by it having no EVs -NEV(!):

In[]:=

Print["EV=",N[ev],", charpolyEV:",cpev,", charpolySV:",cpSigmaSq,", SVsqd=",N[svsqd],", L=",L,", R=",R,", Q=",Q,", Spread=",spread]

EV={0.,0.}, charpolyEV:

, charpolySV:-4x+

, SVsqd={4.,0.}, L=16, R=0, Q=0, Spread=ComplexInfinity

where the “Spread” measure of the matrix elements which is normalized by the linesum EV, ev[[1]], which fails in this rare event where the EVs vanish! Once again the SVs “saved the day”.

RG update!

EDIT update - ResearchGate [52] reported 691 "reads" for my NKS paper [1], and 820 for my paper with Rogers on magic cubes [48]. In the long term I expect those will be eclipsed by the DMPS [4] work with Ian Cameron and Adam Rogers for the information entropy that has been used here in the spectral code.

n=5, rank=5, C=22.7212%, H=1.24375

Acknowledgements

Much of this work began with projects for undergraduate students in my classes - Daniel Schindel, Matthew Rempel, Adam Rogers, Wayne Chan, Marcus Steeds, often later as summer students, and with support from the Winnipeg Foundation. Colleagues Joe Williams, Frank Hruska and Ian Cameron made helpful suggestions. I also had the good fortune to meet and exchange ideas with the late John Hendricks [16] and the late Harvey Heinz [15] (whose site has been kept online), Arno van den Essen [37,41] and Christian Eggermont [36,37], Walter Trump [17], Francis Gaspalou [54], Steve Kirkland [55], as well as recreational mathematicians whom I met online.
Thanks also to George Beck at “The Mathematica Journal” for guidance on Mathematica’s notebook style and the example of Robert Cowen’s article [27] which I used as a guide.

References

[1] (NKS2006) Loly, P.D., “Franklin Squares: A Chapter in the Scientific Studies of Magical Squares”, Complex Systems, 17 (2007) 143-161: https://www.wolframscience.com/conference/2006/presentations/materials/loly-complex_systems-17-1-2.pdf

[2] [IHPST2003] Loly, P.D., “Scientific_Studies_of_Magic_Squares”, International History and Philosophy of Science Teaching (IHPST) Conference, Winnipeg 2003, http://home.cc.umanitoba.ca/~loly/IHPST.pdf

[3] (LAA) Loly, P.D., Ian Cameron, Walter Trump and Daniel Schindel, “Magic square spectra”, Linear Algebra and its Applications, 430 (2009): http://www.sciencedirect.com/science/journal/00243795/430/10 2659-2680. This is an update of the main part of the lead keynote talk at IWMS16, June 1-3, 2007, Windsor, Ontario.

[4] (DMPS) Ian Cameron, Adam Rogers & Peter Loly, “Signatura of magic and Latin integer squares: isentropic clans and indexing”, from 2012 conference published in Discussiones Mathematicae Probability and Statistics, 33(1-2) (2013) 121-149. Download paper: http://www.discuss.wmie.uz.zgora.pl/ps/ - for missing data appendix see next:

[5] Loly, P., “Data Appendix for DMPS” in Mathematica’s bracket notation for DMPS [24]: http://home.cc.umanitoba.ca/~loly/Bedlewo.txt - Wolfram Notebook - DataLoly.nb - Archive c. January 2021 after present Notebook.

[6] Weisstein, E., “Latin Square”, http://mathworld.wolfram.com/LatinSquare.html

[7] Weisstein, E., “Magic Square”, March 3, 2004, http://mathworld.wolfram.com/MagicSquare.html; see also http://mathworld.wolfram.com/notebooks/RecreationalMath/MagicSquare.nb

[8] Swetz, F., “Legacy of the Luoshu: the 4,000 Year Search for the Meaning of the Magic Square of Order Three”, https://turing.cs.hbg.psu.edu/~fjs2/

[9] Onkar Singh, found 14 Jan. 2021 - Contributed by: Onkar Singh (March 2011) -
https : // demonstrations.wolfram.com/MagicSquare/

[10] Nasser, magic[n], https://mathematica.stackexchange.com/questions/73131/the-magic-square-function, “Nasser”, Feb. 3, 2015 (edited March 21, 2018).
https://mathematica.stackexchange.com/users/70/nasser

[11] Cayley, A.C., “Note on magic squares”, Messenger Math. VI (1877), 168.

[12] Dudeney, H., “The magic square of sixteen”, The Queen: The Lady’s Newspaper and Court Chronicle, January 15 (1910) 125-6.

[13] Andrews, W. S., “Magic Squares and Cubes”, The Open Court Publishing Company, 1917, Cosimo 2004.

[14] Benson, W.H. and Jacoby, O., “New Recreations with Magic Squares”, Dover 1976.

[15] Heinz, Harvey, http : // www.magic-squares.net/order4list.htm - see also [9,10] for each of the 48 pandiagonal and associative members at order 4. Also see

[16] Heinz, Harvey D. and Hendricks, John R., "Magic Square Lexicon: Illustrated", 2000 - online at http : // magic - squares.net/Downloads/HendricksBooks/Lexicon - v2.pdf

[17] Trump, W., “Notes on Magic Squares and Cubes”, http://www.trump.de/magic-squares/

[18] RCL) Rogers, A., Cameron, I.D. & Loly, P.D., “Compounding Doubly Affine Matrices”, 2017 arXiv:1711.11084[math.CO]

[19] Weisstein, E., “Panmagic Square”, https://mathworld.wolfram.com/PanmagicSquare.html

[20] See my : http : // home.cc.umanitoba.ca/~loly/PIC1997cover.pdf which has 3, 4, 8 and 20.
Also https : // i.stack.imgur.com/g8wdG.png

[21] Weisstein, E., “Associative Magic Square”, https://mathworld.wolfram.com/AssociativeMagicSquare.html - see also https://en.wikipedia.org/wiki/Associative_magic_square

[22] McClintock, E.(1897) “On the most perfect forms of magic squares, with methods for their production”. American Journal of Mathematics 19 p .99 - 120.

[23] (OllyBrée) Ollerenshaw, K. and Brée, D., “Most-Perfect pandiagonal magic squares: their construction and enumeration”, The Institute of Mathematics and its Applications, Southend-on-Sea, UK, 1998. Note that they used sequences beginning with zero:

0,1,2,..,(

-1),fortheirmagicsquares.

71	64	69	8	1	6	53	46	51
66	68	70	3	5	7	48	50	52
67	72	65	4	9	2	49	54	47
26	19	24	44	37	42	62	55	60
21	23	25	39	41	43	57	59	61
22	27	20	40	45	38	58	63	56
35	28	33	80	73	78	17	10	15
30	32	34	75	77	79	12	14	16
31	36	29	76	81	74	13	18	11

71	8	53	64	1	46	69	6	51
26	44	62	19	37	55	24	42	60
35	80	17	28	73	10	33	78	15
66	3	48	68	5	50	70	7	52
21	39	57	23	41	59	25	43	61
30	75	12	32	77	14	34	79	16
67	4	49	72	9	54	65	2	47
22	40	58	27	45	63	20	38	56
31	76	13	36	81	18	29	74	11

71	64	69	8	1	6	53	46	51
66	68	70	3	5	7	48	50	52
67	72	65	4	9	2	49	54	47
26	19	24	44	37	42	62	55	60
21	23	25	39	41	43	57	59	61
22	27	20	40	45	38	58	63	56
35	28	33	80	73	78	17	10	15
30	32	34	75	77	79	12	14	16
31	36	29	76	81	74	13	18	11

71	8	53	64	1	46	69	6	51
26	44	62	19	37	55	24	42	60
35	80	17	28	73	10	33	78	15
66	3	48	68	5	50	70	7	52
21	39	57	23	41	59	25	43	61
30	75	12	32	77	14	34	79	16
67	4	49	72	9	54	65	2	47
22	40	58	27	45	63	20	38	56
31	76	13	36	81	18	29	74	11

71	64	69	8	1	6	53	46	51
66	68	70	3	5	7	48	50	52
67	72	65	4	9	2	49	54	47
26	19	24	44	37	42	62	55	60
21	23	25	39	41	43	57	59	61
22	27	20	40	45	38	58	63	56
35	28	33	80	73	78	17	10	15
30	32	34	75	77	79	12	14	16
31	36	29	76	81	74	13	18	11

71	8	53	64	1	46	69	6	51
26	44	62	19	37	55	24	42	60
35	80	17	28	73	10	33	78	15
66	3	48	68	5	50	70	7	52
21	39	57	23	41	59	25	43	61
30	75	12	32	77	14	34	79	16
67	4	49	72	9	54	65	2	47
22	40	58	27	45	63	20	38	56
31	76	13	36	81	18	29	74	11